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ABSTRACT 


7  A  Parabolic  Equation  (PE)  workshop  sponsored  by  the  Surveillance 
Environmental  Acoustic  Support  (SEAS)  Project  and  hosted  by 
NORDA's  Numerical  Modeling  Division  was  held  at  NSTL  from  31 
March  to  3  April  1981.  The  purpose  of  the  workshop  was  to  provide 
a  forum  for  those  active  in  theoretical  and  applied  PE  develop¬ 
ment  and  to  compare  computer  results  for  a  set  of  ocean  acoustic 
problems.  Fifteen  formal  presentations  were  given,  and  thirteen 
different  PE  models  were  exercised  against  four  test  cases.  The 
test  cases  ranged  from  one  that  any  PE  should  be  able  to  handle 
to  one  that  is  a  challenge  for  every  PE.  This  report  describes 
the  results  of  that  workshop  and  includes  the  theoretical  bases 
for  the  models,  the  model  descriptions,  results  of  the  test 
cases,  and  abstracts  of  the  formal  presentations.  — 
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NORDA  PARABOLIC  EQUATION  WORKSHOP 


I.  INTRODUCTION 


The  advantages  of  the  parabolic  approximation  (PA)  to  the  elliptic  wave 
equation  for  use  in  obtaining  solutions  to  ocean  acoustic  problems  were  first  made  ' 
known  to  the  underwater  acoustic  scientific  community  through  an  AESD  workshop* 
(sponsored  by  LRAPP)  on  non-ray  tracing  techniques  of  acoustic  propagation  modeling 
held  in  May  1973.  That  workshop  concerned  itself  primarily  with  mode  theory  solu¬ 
tions  to  a  set  of  test  cases  for  range-independent  environments,  but  also  included 
the  Parabolic  Equation  (PE)  results  of  Hardin  and  Tappert^  for  those  same  prob¬ 
lems.  As  is  now  well  known,  the  primary  advantage  of  PE  is  its  applicability  to 
ocean  acoustic  problems  with  range-dependent  environments.  See,  for  example, 
Tappert^  for  the  origin  and  history  of  the  development  and  application  of  PE  to 
various  fields,  especially  the  pioneering  applications  of  Leontovich  and  Fock^. 

Since  that  AESD  workshop  in  1973  there  has  been  a  literal  explosion  of  the 
development  and  application  of  PE  to  solutions  of  underwater  sound  problems.  No 
attempt  is  made  in  this  report  to  try  to  include  what  is  now  an  extensive  PE  bib¬ 
liography.  Because  of  advances  over  the  last  several  years  in  PE  models  for  the 
range-dependent  environment  problems  in  ocean  acoustics,  the  Surveillance  Environ¬ 
mental  Acoustic  Support  Project  (SEAS,  NORDA  Code  520)  sponsored  a  "PE  Workshop" 
hosted  by  NORDA's  Numerical  Modeling  Division  (OSTI.  Code  320).  The  workshop  was 
held  at  NSTL  from  31  March  to  3  April  1981.  Objectives  of  the  workshop  were  to  (a) 
provide  a  forum  for  those  presently  active  in  theoretical  and  applied  PE  development 
to  exchange  ideas,  describe  the  PE  programs,  identify  problems  or  deficiencies  in 
the  PE  approach,  and  to  stimulate  new  ideas  and  approaches;  (b)  compare  computer 
results  for  a  set  of  identical  ocean  acoustic  problems;  and  (c)  collectively  identi¬ 
fy  range  dependent  environment  problems  that  could  provide  a  set  of  benchmark  test 
cases.  An  anticipated  result  of  the  workshop  is  the  identification  of  one  (or 
perhaps  several)  PE  transmission  loss  model  that  could  meet  the  needs  of  the  SEAS 
program.  This  supported  version  would  be  a  "state-of-the-art"  model  derived  from  a 
combination  of  the  most  promising  features  of  existing  PE  models.  The  workshop  was 
a  success  in  tenas  of  achieving  objectives  (a)  and  (b)  above;  however,  objective  (c) 
was  not  achieved — there  was  little  discussion  of  possible  benchmark  test  cases.  This 
report  describes  the  results  of  the  workshop,  but  does  not  include  recommendations 
made  to  SEAS  as  a  consequence  of  it. 

The  workshop  was  divided  into  two  parts:  half  devoted  to  the  test  cases  and 
half  devoted  to  presentations  of  new  i-E  develops  In  the  test  problem  portion  a 

brief  description  of  each  PE  program  exercised  -he  test  cases  was  given  and 

the  predicted  transmission  loss  presented.  Thirteen  different  PE  models  were  exer¬ 
cised  against  the  test  problems,  with  reference  normal-mode  calculations  supplied  by 
the  SACLANT  ASW  Research  Centre.  In  the  second  half  of  the  workshop,  fifteen  half- 
hour  presentations  on  recent  PE  research  and  implementations  were  given.  Attendance 
at  Che  workshop  was  limited  to  invitees.  The  attendees  and  their  organisational 
affiliations  are  given  in  Table  1. 

This  report  addresses  Itself  almost  entirely  to  the  test  case  portion  of  the 
workshop  (Sections  II-IX).  Although  the  new  development  portion  occupied  half  of 
the  workshop  period,  it  was  decided  that  a  semi-detalled  exposition  of  those  presen¬ 
tations  «is  a  major  task  and  beyond  the  scope  of  this  report.  Abstracts  are,  how¬ 
ever,  reproduced  in  Appendix  A. 
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Table  1.  Attendees  and  organizational  affiliation 


NAME 

Dr.  R.  Baer 
Mr.  H.  Brock 
Dr.  H.  Bucker 
Dr.  R.  Cavanagh 
Dr.  S.  Chin-Bing 
Dr.  J.  Davis 
Dr.  L.  Dozier 
CDR.  K.  Evans 
Dr.  R.  Evans 
Dr.  K.  Gilbert 
Mr.  D.  Gordon 
Dr.  R.  Greene 
Dr.  J.  Hanna 
Miss  E.  Holmes 
Dr.  F.  Jensen 
Dr.  D.  Lee 
Dr.  L.  Nghiem-Phu 
Dr.  D.  Palmer 
Prof.  J.  Papadakis 
Mr.  R.  Patton 
Mr.  N.  Paz 
Mr.  M.  Pedersen 
Mr.  J.  Perkins 
Mr.  0.  Seals 
Mr.  C.  Spofford 
Prof.  F.  Tapper t 
Dr.  D,  Thomson 
Mr.  D.  White 
Dr.  0.  Wood 
Dr.  R.  Uenka  tar  email 


ORGANIZATION 


Naval  Research  Laboratory 

Naval  Research  Laboratory 

Naval  Ocean  Systems  Center 

Planning  Systems  Incorporated 

Naval  Ocean  Research  and  Development  Activity 

Naval  Ocean  Research  and  Development  Activity 

Science  Applications  Incorporated 

Naval  Ocean  Research  and  Development  Activity 

Ocean  Dat«  Systems  Incorporated 

Naval  Ocean  Research  and  Development  Activity 

Naval  Ocean  Systems  Center 

Science  Applications  Incorporated 

Science  Applications  Incorporated 

Science  Applications  Incorporated 

SACLANT  ASW  Research  Centre 

Naval  Underwater  Systems  Center 

Daubln  Systems  Corporation 

Ocean  Acoustics  Laboratory  (N0AA/A0ML) 

University  of  Rhode  Island 

Bell  Telephone  Laboratories 

U.  S.  Nival  Oceanographic  Office 

Naval  Ocean  Systems  Center 

Naval  Research  Laboratory 

Bell  Telephone  Laboratories 

Science  Applications  Incorporated 

University  of  Miami 

Defence  Research  Establishment  Pacific  (Canada) 
Naval  Ocean  Research  and  Development  Activity 
Naval  Undersea  Systems  Center 
Bell  Telephone  Laboratories 
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To  assist  In  the  interpretation  of  results  from  the  various  PE's  that  were 
exercised  against  the  test  cases,  theoretical  bases  are  presented  in  Section  II  for 
a  generalized  PE,  standard  PE,  small  modifications  to  PE,  a  C0-independent  PE,  and 
higher-order  corrections  to  PE.  That  section  ends  with  a  summary  of  the  accuracy  of 
the  various  forms.  Section  III  discusses  numerical  algorithms  for  solving  PE,  such 
as  split-step,  finite  difference,  and  finite  element  techniques.  Included  as  well 
are  summaries  of  treatments  of  interface  and  boundary  conditions  and  techniques  for 
obtaining  the  initial  field.  Abstracts  describing  the  PE's  exercised  against  the 
test  cases  are  given  in  Appendix  B. 


Given  the  groundwork  that  facilitates  an  understanding  of  the  component  parts 
of  each  PE  and  allows  proper  perspective  of  their  interrelationships,  the  test  cases 
and  results  are  presented  in  Sections  IV-IX.  Section  IV  discusses  the  rationale  for 
the  test  problems  and  some  "ground  rules"  for  participants.  The  test  cases  them¬ 
selves  are  addressed  in  Sections  V-VIII,  and  each  section  Includes  (1)  description 
of  the  case  and  what  it  was  designed  to  test  and  (2)  model  results,  analysis,  and 
summary.  These  sections,  taken  in  turn,  are  as  follows: 


Section  V. 
Section  VI. 
Section  VII. 
Section  VIII. 


Test  Case  1  -  Range  Dependent  Surface  Duct 
Test  Case  2  -  Bilinear  Profile 

Test  Case  3  -  Range  Independent  Shallow  Water  Environment 
Test  Case  4  -  Basin/Slope/Shelf  Configuration  for  a 
Geoacoustical  Bottom. 


A  summary  of  the  workshop  Is  given  in  Section  IX. 


[  I .  FORMS  OF  THE  PARABOLIC  EQUATIONS 


A.  INTRODUCTION 


Since  there  are  a  number  of  approximations  to  the  Helmholtz  Equation  that  lead 
to  what  are  commonly  referred  to  as  "Parabolic  Equations"  (PE's),  it  is  Important  at 
the  start  to  define  and  distinguish  the  various  forms  presented  at  the  workshop, 
following  Tappert's  “xposlcion  presented  at  the  start  of  the  workshop  and  also 
•  :<>»*cnce  lr  we  introduce  a  generalized  PE  (attributed  to  Claerbout'’*^  and  called 
v'.  .*  here)  from  which  the  various  PE’s  can  be  systematically  derived  and  which 
gives  a  framework  for  identifying  differences  among  them.  See  also  McDaniel^  for 
a  discussion  of  splitting  matrices  and  derivations  of  parabolic  approximations. 

The  treatment  here  will  be  brief,  but  should  help  the  reader  understand  the 
relationships  among  the  various  forms  of  PE.  Such  an  overview  thus  facilitates 
interpretation  of  results  from  the  PE  models  that  were  exercised  against  the  test 
problems. 

In  this  section  then  we  will  derive  a  generalized  PE  and,  in  turn,  derive  stan¬ 
dard  PE,  discuss  small  modifications  to  PE  (CMOD  and  CPA),  outline  the  derivation  of 
a  CQ-iniependent  PE,  derive  higher  order  corrections  to  PE  or  "high  angle”  PE's, 
ami  finally  estimate  the  accuracy  of  the  various  forms. 

B.  GENERAL  PE  (GPE) 

We  begin  by  deriving  a  generalized  PE  following  Tappert-^.  The  acoustic  pres¬ 
sure  p(r,z)  in  a  cylindrical ly  symmetric  medium  of  constant  density  P  having  a  har¬ 
monic  point  source  at  (0,zs)  satisfies  the  following  equation: 
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Prr  +  (Pr/r)  +  Pzz  +  ko2n2P  =  *(Wr)  ps  6(r)6(z  -  zs) 


> 


(1) 


where  r  is  the  radial  variable,  z  the  depth,  Icq  -  w/C0,  oils  the  angular  fre¬ 
quency,  Cq  is  some  arbitrary  reference  sound  speed,  n(r,z)  »  C0/C(r,z)  is  the 
index  of  refraction,  and  p8  the  source  strength.  Eq.  1  is  solved  subject  to 
interface  conditions  of  continuity  of  p  and  p*/p  where  p*  is  the  normal  derivative. 

A  pressure-release  surface  (p(r,0)  »  0)  and  radiation  boundary  conditions  as  z  and  r 
go  to  infinity  are  assumed. 

Away  from  the  source  we  introduce  the  function  u(r,z)  defined  by 

P  -  u/r*5  ,  (2) 


and  make  the  far-field  approximation  kor»l  to  obtain  the  Helmholtz  far-field 
equation 


urr  + 


uzz  +  k02n2u  =  0 


(3) 


Defining  the  operators 


P 


3r 


(4) 


and 

Q  ’  (k?  '  h  * n2)k  <5> 

we  can  then  write 

(P2  +  kG2Q2)u  =>  0  .  (6) 

At  this  point  we  wish  to  “factor”  the  operators  in  Eq.  A.  There  are  various  factor¬ 
izations  (see  Ref.  7  for  example),  but  we  will  use  that  of  Claerbout^  and  write 
Eq.  6  as 


(P  +  ikoQHP  -  ikoQ)u  +  ikolP.Qlu  -  0  , 

where 


[P,Qlu  =  PQu  -  QPu 

is  the  commutator  of  the  operators  P  and  Q.  If  the  medium  is  range  independent,  or 
weakly  range  dependent,  the  commutator  can  be  ignored  and  with  the  assumption  of 
only  “outgoing”  waves  we  arrive  at  GPE: 


Pu  =  ikgQu 


(7) 
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GPE,  according  to  Tappert,  is  the  most  complete  PE  which  is  evolutionary  in 
range  and  neglects  backscattering.  For  range-independent  environments  and  only 
outgoing  waves  it  is  exact  within  the  limits  of  the  far-field  approximation.  It 
provides  the  basis  for  obtaining  the  various  forms  of  "marching  out"  PE’s,  i.e., 
partial  differential  equations  in  first  order  with  respect  to  r  (hence  the  name 
"Parabolic",  as  opposed  to  the  Helmholtz  Elliptic  partial  differential  equation). 
These  forms  will  follow  as  a  result  of  approximations  to  the  pseudo-differential 
operator  Q,  whose  properties  preclude  the  solution  of  GPE  itself. ^  For  later 
convenience  we  also  write  Eq .  7  as 


Pu  =  ikc(l  +  q)^u 


(8) 


with 

Q  -  (1  +  q)*  , 
q  =  e  +  u  , 
e  =  n2  -  1  , 

and 

1  32 

“  ’  k?  ’  a?  • 

It  is  important,  finally,  to  note  that  in  a  range-dependent  environment  GPE 
(and  approximations  to  it)  has  an  intrinsic  error,  [P,Q]u.  Conditions  under  which 
the  commutator  is  small  have  not  been  rigorously  defined,  although  Claerbout^, 
Tappert^,  and  De  Santo^2  provide  frameworks  within  which  co  address  the  problem. 
Until  more  is  learned,  users  of  PE's  must  be  wary  of  cases  in  which  refraction  index 
or  boundary  conditions  are  more  than  weakly  range  dependent. 

C.  STANDARD  PE  (SPE) 

Standard  PE,  which  we  will  refer  to  as  SPE,  is  that  PE  reported  by  Hardin  and 
Tappert  at  the  AESD  workshop2.  Although  It  is  a  first-order,  small  q  approxima¬ 
tion,  SPE  has  the  advantage  of  being  particularly  amenable  to  efficient  numerical 
solution  (viz.,  the  Split-step  algorithm  discuss  below).  It  is  easily  derived  from 
Eq.  8  by  making  a  Taylor  series  expansion  of  Q 


Q  °  (1  +  q2)^  “  1  +  (q/2)  -  (q2/8)  +  (q 3/ 1 6)  +  ••• 


(9) 


and  retaining  only  the  first  two  terras  to  obtain 


Pu  **  ik0(l  +  (q/2) )u  . 


The  Introduc t Ion  of  the  envelope  function 


’{*  (  r ,  z  )  =  u(r,z)exp(-ik0r) 


(10) 
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leads  to 


Pi|>  “  ±(k0q/2)^ 


(ID 


or  to  SPE  in  a  more  familiar  form 

2ik0tpr  +  \pzz  +  k02(n2  -  l)ij>  =  0  .  (12) 

SPE  is  thus  derived  from  the  Helmholtz  far-field  equation  by  assuming  [P,Q]u  =»  0, 
neglecting  backscatter,  and  neglecting  q^  and  higher  order-terms  in  the  expansion 
of  Q  *  vTT  q. 

In  the  usual  derivation  of  SPE  the  envelope  function  is  introduced  in  Eq.  3  to 

yield 

P2ip  +  2ik0Pip  +  k02qij;  =  0  .  (13) 


The  P^^term  is  then  neglected  by  making  the  "parabolic"  approximation,  i.e., 
assuming 


P2ip  |  <<  |2ik0Pij/| 


This  approximation,  which  in  a  special  sense  assumes  that  ^  Is  a  slowly  varying 
function  of  r,  is  implied  by  the  approximations  listed  above  in  obtaining  SPE  from 
GPE.  For  SPE  to  be  a  good  approximation  to  GPE  (in  the  sense  of  some  norm  ||  *  ||) 
the  neglected  terms  in  the  expansion  of  (1  +  q)*5  mast  be  small,  i.e., 

IMF  »  IF  Ml2 

must  be  small  compared  to  one. 

Modal  analysis  by  McDaniel^  for  a  range-independent 
given  the  correct  field  at  some  range  rQ,  the  mode  depth 
identical  to  those  for  the  far-field  Helmholtz  equation, 
range  r  for  the  n-th  mode  depends  on  ICq  and  is  given  by 

((kQ2  +  kn?)/2ko)(r  -  r0)  f 

as  opposed  to  k^r  -  rQ)  for  the  far-field  Helmholtz  equation.  The  SPE  phase  can 
.1I90  be  wrttten  in  terms  of  an  equivalent  angle  9n  with  respect  to  the  horizon¬ 
tal  as 

(k0/2)(l  +  cos20n)(r  -  rD)  , 


environment  shows  that 
functions  for  SPE  are 
However ,  the  phase  at 
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where 


kn  =  ^cosc^  . 


The  ratio  of  the  SPE  mode  function  to  the  correct  mode  function  for  the  n-th  mode  is 
therefore  given  by 


RSpe  =  exp(i2k0s2 (r  -  r0)) 


(14) 


where 


s  =  (cin(0n/2) ) 2 


(15) 


Motice  that  the  phase  error  in  Eq.  14  depends  on  the  choice  of  k0.  One  can 
select  k0  such  that  the  phase  err<v:  is  "oro  for  a  sing!  »  irode,  t.e.,  k0  »  kjj, 
hut  not  for  all  modes.  If  the  modes  are  contained  in  a  sufficiently  small  band 
(i.e.,  the  spectrum  of  the  operator  [kgQ]  is  concentrated  near  kQ),  k0  can  be 
selected  to  be  in  the  center  of  the  hand  and  the  resulting  phase  errors  and  their 
consequences  on  the  field  may  be  small.  This  so-called  "narrow  band"  case  is  usual¬ 
ly  interpreted  as  a  "small  angle"  restriction;  the  horizontal  components  of  the  true 
wave  numbers  kjj  »  k0cos8n  are  confined  to  a  small  band  about  l^,  and  the 
equivalent  geometric  rays  are  nearly  horizontal  (8n  as  0).  As  noted  here  and  by 
Fitzgerald^  and  Brock,  Buchal  and  Spoffard*^  the  angles  need  not  be  small,  but 
rather  limited  to  a  small  aperture  so  that  k^cosQj,  is  nearly  constant.  For  a 
wide  band  the  phase  errors  can  be  large  with  significant  consequences  on  the  field. 

The  next  two  ?E  forms  discussed  (CMOD  and  CPA)  are  intended  to  reduce  phase 
errors  while  retaining  those  features  of  SPE  which  allow  for  efficient  computer 
solution. 

0.  CMOD 

A  modification  to  SPE  that  is  designed  to  reduce  the  associated  phase  errors  or 
improve  the  narrow  band  limitations  has  been  developed  by  Brock,  Buchal  and  Spof- 
ford  •  This  modification,  referred  to  as  CMOD,  Is  not,  insofar  as  we  know,  easily 
derivable  from  GPE.  SPE  phase  errors  have  been  viewed  as  the  result  of  errors  in 
each  mode's  horizontal  wave  number  (kn),  which  can  be  thought  of  as  errors  in  each 
mode's  phase  velocity  (w/kp).  CMOD  is  based  on  the  idea  of  constructing  a  “pseudo 
problem  whose  (SPE)  phase  velocities  are  equal  to  the  elliptic  phase  velocities  of 
the  corresponding  modes  in  the  original  problem"*'1.  This  is  accomplished  by  util¬ 
izing  the  UKB  approximation  to  identify  a  mapping  of  the  modal  turning  point  depths 
tn  the  original  problem  into  the  turning  po  int  depths  of  the  same  mode  in  the  pseudo 
problem.  Carried  out  exactly  this  produces  a  zero  phase  error,  but  rhe  depth  func¬ 
tions  are  then  in  error  away  from  the  mode  turning  points.  If  the  main  contribution 
to  the  field  of  any  mode  is  at  a  depth  in  the  vicinity  of  its  turning  point,  then 
this  will  be  a  small  e-ror  and  will  be  independent  of  range  in  a  range- independent 
environment,  bse  of  the  WKB  approximation,  however,  does  not  reduce  the  phase  error 
to  exactly  zero;  therefore,  there  are  small  phase  errors  as  well  as  depth  function 
errors . 
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(16) 


CMOD  employs  a  frequency  independent  mapping,  given  by 
(n,z)  (m,>)  =  ((2n  -  1)^,  z\-H)  . 

With  this  mapping  the  form  of  the  PE  equation  is  unchanged,  i.e.,  u(r,y)  satisfies 


Pu  =  il.c[ 


m2  +  1 


1  8' 


2  2k02  3y2 

or  in  terms  of  the  envelope  function 


•]u  , 


(17) 


2ik0>f'r  +  4>yy  +  k02(m2  -1)^=0  , 


(18) 


which  is  the  same  form  as  SPE  given  by  Eq.  12. 

Analysis  of  phase  and  amplitude  errors  of  SPE  using  CMOD  is  complicated  even 
for  the  range-independent  environment  case.  See  Reference  3.0  for  details  and  cau¬ 
tions  regarding  its  use  in  problems  having  isovelocity  regions. 

E.  CPA 

Palmerll  showed  that,  in  a  range-independent  environment,  an  integral  opera¬ 
tor  applied  to  parameterized  solutions  of  SPE  yielded  a  solution  to  the  Helmholtz 
equation.  The  SPE  solution  was,  in  fact,  found  to  be  the  result  of  a  stationary 
phase  approximation  to  the  integral.  An  expansion  about  the  stationary  phase  point 
yields  a  relation  of  form 


4>(r,z)  =  {1  -  (q/4)  +  (q2/8)  +  ••• 

+  ik0r  (-(q2/8)  +  •••) 

+  (1/12) (ik0r)2(q/2)4  +  — }*(r,z)  ,  (19) 


where  $  is  the  SPE  envelope  and  <t>  is  the  Helmholtz  equation  envelope.  Note  that  the 
first  and  fourth  terms  listed  give 

<J>(r,z)  K  (1  -  i(k0rq2/8))iKr,z) 

(20) 

-  i|>  +  (ir/2k0)'l'rr  . 

De San  to used  a  similar,  but  different,  integral  operator  to  relate  4>  and  ^  for  a 
general,  range-dependent  environment.  Again,  expansion  about  the  stationary  phase 
point  leads  to  a  result  much  like  Eq.  19,  but  includes  the  effects  of  range-depend¬ 
ent  n(r,z) .  DeSanto,  Perkins,  and  Baer  suggested  using  Eq.  20  to  improve  SPE,  and 
showed  examples  of  the  improvement This  modified  solution  is  termed  the  Cor¬ 
rected  Parabolic  Approximation,  or  CPA,  and  takes  the  form,  At  each  range  step, 


ii>„(r)  *  iKr)  +  (iA/2KQHr>.(r) 


(21) 
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where  ^(r)  is  the  SPG  solution,  A  is  the  range  step,  and  v'rr(r)  is  a  finite- 
difference  approximation  to  ^rr  which  utilizes  ^(r-A),  r),  ^(r+A).  See  Ref.  13 

for  details  of  the  implementation. 

CPA  can  be  derived  directly  from  GPE  via  the  second-order  expansion  of  Q: 


p* 2  ik0((q/2>  -  (q2/sm 


(22) 


and  the  approximations : 


^rr  “  “(k02q2/^)^» 


for  q  k  0. 


(23) 


(r  +  A)  *  e 


ik„((q/2)  -  (q2/8))A 


•  ip(r) 


(24) 


for  q,  q2,  and  qr  .small.  Tappert  emphasized  at  the  workshop  that  this  modifica¬ 
tion  is  valid  only  for  small  q,  in  which  case  the  correction  is  small  anyway.  A 
phase  error  analysis  below  confirms  and  quantifies  this.  He  suggested  that  a  more 
useful  (more  accurate)  solution  to  second  order  in  q  is  obtained  via  Claerbout's 
rational-linear  method  (see  Subsection  G). 

Modal  analysis  of  CPA  for  a  range-independent  environment  leads  to 

Repa  *  (1  +  Y2)^exp  i(2k0s2(r  -  r0)  -  tan-1y)  >  /2 


where  R  is  the  ratio  of  the  CPA  mode  function  to  the  Helmholtz  mode  function. 


y  «  2k0(r  -  r0)(s  -  s2)2 


(26) 


and  as  before  s  ■  sin2(8n/2).  Thus,  for  the  general  case  CPA  contains  an 
amplitude  and  a  phase  error.  Notice  that  the  first  term  in  the  exponential  is  the 
SPE  phase  error,  Eq.  14,  and  the  second  term  is  a  modification  to  that  phase  error. 

For  small  s  (small  8n,  and  only  the  lowest-order  terms  in  amplitude  and 
phase) ,  Rcpa  reduces  to 

R cpa  %  +  2ko2s4^r  ~  r0)2)exp(i.41t0s3(r  -  rQ))  .  (27) 

s-+o 

In  this  event  CPA  provides  a  more  accurate  phase  than  SPE  at  the  expense  of  a  small 
amplitude  error.  Note,  however,  that  when  s  is  not  small  (i.e.,  the  "wide  band” 
case)  the  error  is  not  small.  This  is  the  case  for  any  small  q  approximation. 
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F.  C0-INDE PENDENT  PE 

As  has  been  noted  in  Section  C,  SPE's  accuracy  depends  on  | |e  +  m|  |  being 
small.  Tappert^  has  shown  however  th3t  it  is  possible  to  improve  on  SPE  by  either 
(i)  requiring  ||e||  to  be  small  but  relaxing  the  condition  on  n  and  thereby  allowing 
wide  beamwidths,  or  (2)  requiring  ||/»||  to  be  small  but  relaxing  the  condition  on  e 
and  thereby  allowing  large  variations  in  the  index  of  refraction.  Tappert  uses 
formal  operator  expansions  for  Q  that  are  linear  in  one  of  u.  or  e  and  valid  to  all 
orders  in  the  other.  His  result  for  small  ||/i||  is 


1 


nz' 


Pu  =  iko  [n  +  3 

4k0  n3 


-  +  _L_ 

nz  2k02  n  oz 


(28) 


This  is  Tappert’ 8  CQ-  or  ^-independent  PE  equation.  That  it  is  Icq -independent 
is  easily  seen  since  kQ  and  n  always  appear  in  the  combination  kQn. 

Modal  solution  of  Eq.  28  for  a  range-independent  environment,  to  compare  with 
the  modal  solution  of  the  far-field  Helmholtz  equation,  is  complicated  by  the  fact 
that  the  depth  functions  are  different,  and  we  shall  not  attempt  that  analysis  here. 
If,  however,  the  z  derivative  of  n  can  be  ignored,  then  Eq.  28  reduces  to  SPE  and 
the  corresponding  SPE  phase,  errors  persist. 

Tappert  also  derives  a  small  j|ej|  PE,  but  since  it  was  not  exercised  against 
the  workshop  test  problems  we  refer  the  reader  to  Reference  3. 

G.  RATIONAL  LINEAR  PE  APPROXIMATION 

In  contrast  to  Tappert's  small  ||m||  or  small  ||e||  improvement  over  the  linear 
q  representation  of  the  operator  Q  that  results  in  SPE,  there  are  higher  order 
approximations  of  Q  that  are  correct  to  order  q2.  These  are  generally  referred  to 
as  "Rational  Linear"  approximations,  and  take  the  following  form 


Q  =  (1  +  q)1^ 


a  +  bq 
c  +  dq 


(29) 


To  obtain  SPE  from  this  form,  let  a=l,  b=l/2,  c=l,  and  d=0.  Claerbout^’^,  who 
apparently  was  the  first  to  use  a  rational  linear  approximation  for  PE,  has 
suggested  the  values:  a=l,  b«3/4,  c=l,  and  d=l/4,  so  that  Q  is  approximated  as 

0  =  Jl_±  „(3q/4) 

W  1  +  (q/4) 

Expansion  of  the  denominator  of  Eq.  30  leads  to 


Q  =  1  +  (q/2)  -  (q2/8)  +  (q 3 / 32 )  +  ••• 


(31) 


and  comparison  of  Eq.  31  with  the  Taylor  series  expansion  of  Q  given  by  Eq.  9  shows 
that  it  is  correct  to  order  q2  and  is  therefore  t.  higher-  ordet  (also  higher- 
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angle)  approximation  to  Q  than  employed  to  obtain  SPE.  Use  of  Q  given  by  Eq.  30  in 
5q.  7  leads  to  the  governing  equation  for  u: 


(1  (q/4))Pu  =  ikQ(l  +  (3q/4))u  . 


(32) 


Notice  that  this  is  linear  in  q,  i.e.,  a  partial  differential  equation  that  is  of 
second  order  in  the  z  derivative,  whereas  retention  of  the  q2  term  in  the  Taylor 
series  expansion  of  Q  given  by  Eq.  9  would  result  in  an  equation  of  fourth  order  in 
the  z  derivative.  If  we  introduce  the  envelope  function  Claerbout's  rational 
1  inea-:  approximation  becomes 


(q/4)<J>r  +  -<|>r  =  (ikQq/2)4i  .  (33) 

If  the  first  term  in  Eq.  33  is  set  equal  to  zero,  the  result  is  SPE,  Eq.  12. 

Claerbout's  PE,  Eq.  33,  can  be  arrived  at  more  directly  from  Eq.  13  by:  (1) 
initially  neglecting  the  P2  ¥  term,  (2)  applying  P  and  solving  for  P2  ^  ,  and 
(3;  back  substitu  *ng  the  P2  ^  term  in  Eq.  13.  This  rationale  is  employed  by 
Thomson  (this  report)  and  is  given  in  Ref  &. 

Greene  (this  report)  has  suggested  another  rational  linear  approximation  to  Q 
given  by 


U. 99986  +  0. 79624q 
1.0  +  C . 309Sq 


(34) 


where  the  numerical  coefficients  (a,h,e,  and  d)  have  been  determined  by  minimizing 
the  maximum  error  of  the  phase  of  the  propagating  modes  over  a  range  of  angles  from 
0  to  40°,  i.e.,  kn/k0  between  1.0  and  0.77.  Notice  that  Greene's  coefficients 
in  Eq.  34  are  nearly  the  same  as  Claerbout's  <n  Eq .  30.  Application  of  Greene’s 
mini-max  technique  for  a  range  of  angles  other  than  0-40°  will  of  course  alter  the 
numerical  values  of  a,  b,  c,  and  d. 

Any  rational  linear  approximation  (i.e.  of  form  given  by  Eq,  29)  leads  to  the 
following  envelope  function: 


o^r  +  ’Y'l'zzr  =  +  mij/zz  ,  (35) 

where 

u  =  c  +  d(n2  -  1)  , 

8  =  ik0((a  -  c)  +  (b  -  d)(q2  -  1))  , 
y  =  d/k02  , 
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and 


to  =  i(b  -  d)/k0 

Eq.  35  avoids  the  numerical  specifics  and  introduces  a  convenient  notation  for  a 
later  section  of  this  report. 

Note  that  the  rational  linear  approximation,  regardless  of  the  choice  of  a,b,c, 
and  d,  has  the  same  mode  depth  functions  as  the  far-field  Helmholtz  equation  in  a 
range-independent  environment  (as  would  any  order  Taylor  series  expansion  of  Q),  and 
yet  is  second  order  in  the  z  derivative  in  range-dependent  environments  (unlike  a 
higher-order  Taylor  series  expansion  of  Q). 

Modal  analysis  for  a  range-independent  environment  for  Claerbout's  PE,  Eq.  33, 
results  in  the  error  ratio: 


Rcb  =  exp(i2k0(r  -  rc)s3/(l  -  s  +  s2)) 


(36) 


where  again  s  =  sin2(0/2)  and  kjj  =  k0cosQn.  Notice  that  the  amplitude  is  correct 
and  that  the  phase  e’Tor  is  of  order  as  s-K),  while  the  SPE  error  is  of  order 
s^  (Eq.  14).  The  error  in  Greene's  approximation  is  of  order  1  as  s-K),  but  is 
consistently  "small"  over  the  range  of  angles  0  ^  0n  _<  40° . 

H.  SUMMARY  OF  PE  FORMS  AND  ACCURACY 


This  section  (II)  has  endeavored  to  logically  relate  various  PE  forms  to  a  gen¬ 
eral  PE,  which  in  turn  is  closely  related  to  the  Helmholtz  equation.  GPE  (in  the 
far-field  approximation)  is  exact  when  the  environment  is  range-independent  and  only 
"one-way"  propagation  is  present.  If  nr  ^  0  then  errors  are  made.  Tappert^ 
gives  a  qualitative  ranking  of  those  errors,  and  DeSanto^  implicitly  estimates 
them  in  relating  SPE  solutions  to  solutions  of  the  Helmholtz  equation,  but  rigorous 
quantitative  estimates  remain  to  be  derived.  The  error  could  be  large  whenever  a  PE 
is  used  and  must  be  kept  in  mind. 


The  next  type  of  error  is  that  incurred  in  making  approximations  to  GPE.  GPE 
is  not  amenable  to  numerical  solution;  the  operator  (1  +  q)^/^  is  "non-local"  and 
cannot  be  represented  by  a  finite  sum  of  linear  differential  operators.  The  most 
popular  approximation  to  GPE  to  date  has  been  "SPE",  a  small  q  expansion  and,  as  we 
discuss  in  Section  III,  its  widespread  use  is  a  result  of  an  efficient  numerical 
algorithm  (Tappert  and  Hardin's  split-step)  for  its  solution. 


Error  estimates  of  SPE  and  other  PE's  have,  to  date,  concentrated  on  modal 
analysis  for  range-independent  environments.  Although  such  an  approach  ignores  the 
effect  of  environmental  range  dependencies,  it  is  useful  in  estimating  local  errors 
and  their  cumulative  effect.  As  outlined  in  Section  C,  SPE  has  phase  errors. 
Because  of  this,  a  number  of  PE's  have  been  developed  to  improve  on  SPE. 


(1)  CMOD  -  A  physically  intuitive  correction  to  SPE  whose  quantitative 
improvement  over  SPE  is  not  easy  to  calculate,  but  which  has  been  demonstrated  in 
some  deep  water  ocean  acoustic  problems^. 
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(2)  CPA  -  Error  estimates  here  indicate  an  improvement  over  SPE  for  small 
angle  propagation.  It  has  a  demonstrated  advantage  for  some  realistic  problems^- 3. 

(3)  Tappert's  C0-Independent  PE  -  An  improvement  over  SPE,  but  not  amenable 

to  error  analysis.  Tests  have  not  been  documented  (other  than  presented  here  in  the 
test  problem  sections),  and  the  method  may  have  limited  applicability  since  large 
variation  in  the  index  of  refraction  (large  e  )  is  usually  important  for  wide  angle 
propagation  (large  fi). 

(4)  Rational  Linear  Approximation  -  Both  the  Claerbout  and  Greene  approxi¬ 
mations  offer  the  potential  for  a  significant  improvement  over  SPE — for  "wide  angle” 
propagation.  In  modal  terms,  the  correct  depth  functions  are  retained,  while  the 
phase  errors  are  significantly  reduced  for  angles  up  to  40°  or  more. 

m.  METHOD  OF  SOLUTION 

A.  INTRODUCTION 

The  primary  numerical  advantage  of  the  parabolic  approximation  is  that  it 
yields  a  first-order  differential  equation  in  the  range  variable  (r),  and  can 
therefore  be  solved  by  marching  out  in  range  from  a  given  "initial"  field.  This 
section  is  then  devoted  to  a  discussion  of  numerical  methods — techniques  for  solving 
a  given  PE  form  on  a  computer.  In  practice  the  solution  is  obtained  by  laying  down 
a  numerical  grid  in  depth  z,  solving  for  the  field  at  the  grid  points,  and  advancing 
in  range  to  the  next  range  increment.  The  two  basic  types  of  approaches  are  the 
split-step  algorithm  and  classical  finite  difference/finite  element  techniques.  An 
integral  part  of  each  approach  is  its  treatment  of  interfaces/boundaries,  approxima¬ 
tions  to  the  radiation  conditions,  and,  of  course,  the  initial  field. 

In  this  section  we  address  Split-Step  SPE  (Part  B),  finite  difference  methods 
for  SPE  (Part  C) ,  and  a  comparison  of  the  two  methods  (Part  D).  Part  E  describes 
approaches  for  high  angle  (rational  linear)  PE  and  finally  Part  F  presents  various 
ways  to  obtain  the  initial  field. 

B.  SPLIT  STEP  SPE 

The  split-step  algorithm  was  historically  the  first  technique  employed  in 
solving  SPE  for  underwater  acoustic  problems;  it  is  therefore  appropriate  to  start 
with  this  approach.  To  facilitate  comparisons  among  numerical  solutions  to  SPE,  we 
first  introduce  the  operators: 


x  =  (ik0q/2)  =  A  +  B 


07) 


with 

■  i  32 
=  2k0  3z2 

and 


08) 


B  =  (ik0/2) (n2  -  1) 


(39) 


13 


Write  SPE,  Eq.  11,  as 


=  (A  +  B)iJj 


(40) 


where  the  prime  is  understood  to  be  the  r  derivative. 

Now,  given  the  value  of  ^  at  some  range  rn,  i.e.,^n  =*  ^(r n,z),  the  value 
of  f  at  a  new  range  r+A,  i.e.,iAn+i  »  ^(rn+A,  z),  can  be  approximated  from  Eq.  40 
as 


xA 

e 


^  +  E 

n 


(41) 


where  x  approximates  the  behavior  of  x  over  the  interval  A,  and  where  E  is  an  error 
term  (equal  to  zero  if  X^x  and  x  is  r  independent).  The  actual  error  depends,  as  we 
shall  see  later,  on  how  the  exponential  operator  is  split  and  on  how  x  is  approxi¬ 
mated  by  x  in  a  range-dependent  environment. 

The  split-step  algorithm  usually  takes  one  of  two  forms.  The  first  approach, 
originally  proposed  by  Tappert^,  is  to  split  the  exponential  operator  as  follows: 


e(A  +  B)A  ~  eBA  .  eAA 


(42) 


and  exploit  the  Fourier  transform  to  calculate  the  two  z  derivatives  inherent  in 
exp  (AA) .  is  thus  found  by  means  of  the  algorithm 


,BA 


_  is2A 
2k0 

(e  F  (<J;n) ) 


(43) 


where  F  is  the  Fourier  transform  from  z  to  8  and  F-*  is  the  inverse  transform  from 
8  to  z.  The  second  approach,  later  suggested  by  Tappert^,  ia  to  split  the 
exponential  operator  as 


e(A  +  B)A 


AA 
e  2 


BA 


AA 
e  2 


(44) 


and  solve  for  by  the  algorithm 

1 2A  _^s2A 

i  —  r- 1  r  4k0  *. ,  BAr  —  i  /  4ko^-  ,  . .  ,  (**3) 

♦n+l  ■  F  Me  F{e  F  Me  Ff^))})  , 

In  either  case  appropriate  boundary  conditions  are  required  in  order  to  numeri¬ 
cally  implement  this  approach.  The  physical  conditions  are  that  p-0,  or  ^=0,  at  z»0 
and  ^-*0  as  Since  one  cannot  numerically  deal  with  an  infinite  transform,  the 

latter  (radiation)  boundary  condition  must  be  approximate!.  The  usual  approach^ 
is  to  add  an  attenuating  layer  of  thickness  D  to  the  physical  depth  H  and  set. 
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at  .  The  attenuating  layer  Is  modeled  with  a  complex  index  of 

refraction  of  the  form 


n2(z)  =  n2(H)  +  ia(z);  H  <  z  <  z 

—  max 


where  «(z)  is  the  assumed  attenuation  profile  in  the  layer.  The  choice  of  D  and  a  is 
arbitrary  only  to  the  extent  that  no  significant  energy  should  return  into  the  prob¬ 
lem  from  this  "false"  bottom.  Usually  D  is  taken  to  be  H/314,  i.e.,  zmax  -  4H/3, 
and  a  is  exponentially  Increasing  with  depth.  th  this  approximation  to  the  radia¬ 
tion  condition,  the  Fourier  transforms  in  Eqs.  43  and  45  can  be  replaced  by  a  dis¬ 
crete  FFT  (in  fact,  hy  a  fast  sine  transform  since  p=0  at  z=0  and  2=zmax;).  The 
advantages  of  the  split-step  FFT  are  that  it  is  energy  conserving,  unconditionally 
stable,  and  computationally  efficient2.  The  disadvantages  are  that  It  requires  a 
uniform  z-grid,  periodic  boundary  conditions  and,  as  we  shall  see,  small  step-size 
for  large  propagation  angles  (subject,  of  course,  to  the  limits  of  SPE  itself).  In 
order  to  handle  realistic  ocean  acoustic  problems  rapid  variations  ("discontinui¬ 
ties  )  in  sound  speed,  attenuation  and  density  must  be  taken  into  account.  Sound 
speed  anil  attenuation  are  both  included  in  the  index  of  refraction,  and  the  remedy 
is  simple  but  costly:  model  discontinuities  by  a  continuous  change  over  a  relative¬ 
ly  small  depth  interval  and  force  a  larger  transform  size  for  proper  sampling  In 
depth. 


Density  variations  or  discontinuities  are  a  special  problem.  Tappert3  has 
shown  that  if  one  introduces  a  new  variable  into  Eq.  1 

**  —  k 

P  “  PP  ^ 


and  proceeds  as  in  Section  C  with  u=*p  r3/2,  exp(-ik0r),  the  resulting 
equation  for  f  is  exactly  of  the  form  of  SPE  for  ^  (Eq.'  11),  except  that  n2  (which 
can  also  contain  absorption)  is  replaced  by  n2  where 

fi2  51  °2  +  2k^  *  ^  *  7  *  (P*3/2VP)  . 

Large  gradients  or  discontinuities  in  p  can  make  n  deviate  considerably  from  unity 
(thus  violating  conditions  for  validity  of  SPS).  Tappert  concludes  that  a  “smearing 
out  of  such  rapid  or  discontinuous  changes  is  necessary  and  suggests  an  analytic 
form  for  succ  a  smearing3.  Quite  apart  from  the  particular  analytic  form,  the  end 
result  Is  a  larger  transform  size  and  smaller  range  step. 


C.  FINITE  DIFFERENCE  FOR  SPE 


An  alternative  nunerical  technique  for  solution  of  SPE  Is  provided  by  finite 
difference  techniques.  A  primary  advantage  of  this  approach  is  that  it  does  not 
require  a  uniform  z  mesh  and  discontinuities  in  sound  speed  and  density  can  be  more 
accurately  approximated  without  exorbitant  Increases  in  mesh  points. 

A  standard  finite  difference  approach,  Crank-Nicholson,  applied  to  SPE  (Eq.  II) 
for  constant  z  and  r  rcosh  Has  the  following  forn; 
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(47) 


n 


n  .  n 


1  n+1  n  1  n+1  n+1  n+1  „  „  „  v , 

T  (<h  '  4'n  )  +  TTjt  (V  -  2<l>  +  *P  )  +  (i[>  -  2i|>  +  )  ] 

A  z  1  2  62  £+1  A  £-i  £+1  A  A-l 


k02r/  0  ^+1  11+1  /  0  ,  n, 

+  -?-[(n2  -  1)  +  (n2  -  1)  1|»  ]  =  0  , 

4  £  £  £  £ 


where  £  is  the  z  mesh  index  with  interval  8,  and  n  is  the  r-mesh  index  with  interval 
A.  By  rearranging  Eq.  47  and  introducing  $  as  the  column  vector  whose  components 
are  then's,  we  obtain  the  matrix  form 


—  n+1 

°L  ♦  =  CR 


_  n 


(48) 


Since  Cl  and  Cg  are  tridiagonal  matrices  (even  for  variable  z  and  r  mesh),  Eq. 

48  is  easily  solved  implicitly. 

Crank-Nicholson  is  just  one  of  many  finite-difference  schemes  applicable  to 
SPE.  Others  offer  different  orders  of  accuracy,  ease  of  implementation,  and  speed. 
In  any  of  these,  attenuation  can  be  incorporated  by  use  of  a  complex  index  of 
refraction.  Rapid  variations  in  sound  speed  and  density  are  easily  accommodated  by 
the  variable  mesh  (8)  and  the  flexibility  to  satisfy  the  interface  conditions,  i.e. 
continuity  of  pressure  p  and  p'/p  ( p *  is  the  normal  derivative). 

Radiation  conditions  as  z->-  are  typically  handled  as  for  split~step,  i.o., 
with  a  false  bottom.  An  exception  to  this  is  the  use  of  impedance  boundary 
conditions  by  Papadakis  (this  report),  given  in  terms  of  an  integral  along  the 
bottom  interface 

r  kQ  „ 

.*5  i-r-  (n  2  -  1)  (r  -  s) 

i  PB  /  2  B  L 

p(r,zfi)  =  -  - - *  /  e  *(r  -  s>^  Pz(s,zB)ds  , 

°"P  J  (49) 

o 

where  zg  is  the  interface  depth  below  which  there  is  a  half  space  with  density 
Pg  and  index  of  refraction  ng.  This  technique  avoids  the  introduction  of  a 
non-physical  subbottom  layer  and  offers  a  possible  savings  through  elimination  of 
about  25X  of  the  mesh  points  and  evaluating  Eq.  49  instead. 

0.  SPLIT-STEP  VS.  FINITE  DIFFERENCE  FOR  SPE 

In  this  section  we  will  compare  the  range-step  dependent  and  frequency 
dependent  errors  incurred  in  using  the  split-step  algorithm  and  finite  differencing 
to  solve  SPE.  Included  In  the  results  are  estitaates  made  by  Brock*4  Buchal  and 
Tappert^  McDaniel2,  DiNapoll  and  Deavenport ^2,  and  Jensen  and  Krol^®. 

Assume  chat  the  value  of  *  ^(r„,z)  is  available  at  range  rn.  A  power 
series  solution  to  SPE  in  the  fora  of  Eq.  40  is  given  to  third  order  in  A  by 

<i>  ,  »  U  +  xA  +  <x'  +  x2)|2  +  (x"  +  2x'x  +  XX*  +  x3)|3]  *  ,  (50) 

n+1  2  o  n  n 
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where  the  prime  indicates  an  r  derivative  and  x,  x'  and  x"  are  evaluated  at  rn. 

Consider  first  the  errors  inherent  in  the  exponential  solution  of  SPE  defined 
by  Eq.  41,  i.e.,  even  if  x“x,  the  solution  is  in  error  whenever  xVO.  If  x  is  x 
evaluated  at  rn,  then  an  expansion  of  the  exponential  for  small  A  leads  to  the 
local  error  (difference  between  the  solution  of  Eq.  44  and  Eq.  <-.5): 


2 

E(x(r  ))  =  4  x  '  q,-  +  0(A3)iJ/ 

u  L  n  n  r 


(51) 


whereas  if  x  is  x  evaluated  at  rn+Y>  the  error  is  (to  third  order  in  A): 


E(x(rn  +  p)  =2 3  (x"  -  2xx*  +  Zx'x)^ 


(52) 


Thus,  errors  of  order  A^  can  be  incurred  even  before  the  exponential  operator 
exp  (A+B)  is  split.  These  errors  persist  in  trtiat  follows. 

Now  consider  the  split-step  algorithm  given  by  Eq.  42.  If  x  is  x  evaluated  at 
rn  the  error  is 


Eg’B(x(rn))  =  j  (x'  +  [A.B])^  .  (53) 

This  is  of  order  A^  and  depends  on  the  range  dependence  of  x  (through  n)  and  on 
the  commutator  of  A  and  B.  Comparing  Eqs .  51  and  53  we  see  that  the  splitting 
introduces  an  additional  error,  involving  [A,B].  Excluding  boundary  considerations 

[A,B]*n  =-{(n2)zz*n  +  2(n2)z(«i,n)zz  .  (54) 


The  magnitude  of  [A,B]  is  therefore  seen  to  depend  on  the  z  gradient  of  the  index  of 
refraction  and  is  independent  of  frequency.  On  the  other  hand 


x 


» 


» 


and  is  linear  in  frequency 
to 


A 

x  evaluated  at  rn  +  ^  with  the  result 


The  x'A^  error  term  can  be  removed  by  taking  5  equal 


,,A,B 

Lss 


(x(r 


“  (A,B) 


n*n 


(55) 


This  is  to  be  compared  with  Eq.  52  where  the  operator  was  not  split. 

Consider  now  the  splitting  given  by  Eq.  44  with  x  equal  to  x  evaluated  at  rn. 
The  resulting  error  is 


(*<rn)) 


» 
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(56) 


identical  to  Eq.  51  and  an  improvement  over  the  previous  splitting  but  still  of 

order  A2.  If  on  the  other  hand  x  is  x  evaluated  at  rn  +  ^  the  error  becomes 

“  2 


E *2(x(r  +  A))  =  A-  (B"  +  2(B*A  -  AB')  +  [A, [A,B] ]  +  2[B, [A,B] ]) J  (57) 
ss  n  2  24  n  n 


and  is  now  of  order  A^.  Notice  that  the  term  containing  the  commutators  is  explic¬ 
itly  linear  in  frequency  as  are  other  terms.  Also  notice  that  since  A  is  not  a  func¬ 
tion  of  r,  Eq.  57  differs  from  Eq.  52  only  in  the  commutator  terms.  These  additional 
error  terms  are  a  consequence  of  splitting  the  operator  exp  (xA). 

Before  proceeding  further  it  is  important  to  understand  that  the  error  in  the 
above  formulations  and  those  we  will  see  later  depend  on  the  r  dependence  of  n.  If 
n'  is  small  the  error  will  be  small,  and  if  n'  is  large  the  error  will  be  large — 
thereby  requiring  a  small  range  step.  Regardless  of  the  size  of  n*  one  can  always 
choose  A  to  make  the  local  error  as  small  as  one  likes.  However,  SPE  (as  well  as 
GPE)  is  based  on  the  assumption  that  n*  is  small,  i.e.,  [P,Q]  *  0.  Thus,  in  a 
problem  where  n’  drives  the  range  step-size  error,  SPE  may  be  a  poor  approximation. 

We  next  consider  the  finite  difference  approaches.  One  scheme  for  solving  Eq. 
40  is  given  by 

tl;  il)  iL 

n+1  ~  n  =  n+1  +  n  (58) 

A  *n  2  ’ 

and  results  in  an  error 


EFD 


x'A2 

2 


(59) 


i.e.,  of  order  A2,  linear  in  frequency,  and  identical  to  Eqs.  51  and  56.  However, 
use  of  the  Crank-Nicholson  algorithm 


K+l  -  K  *n+l  »n+l  +  Vn 
A  "  °  2 


(60) 


results  in  an  error  given  by 


E 


CN 


A3 

— -  (x"  +  2x'x  +  xx* 
12 


(61) 


and  is  of  order  A^  and  cubic  (through  x^)  in  frequency. 

We  have  seen  that  the  accuracy  of  split-step  or  finite  difference  algorithms 
depends  ou  how  the  range  dependence  of  x  is  approximated.  For  x  simply  evaluated  at 
range  steps,  the  best  that  one  can  do  is  order  a2,  whereas  a  linear  approximation 
for  x  between  range  points  allowr  order  A^.  Therefore,  in  summary,  in  a  range- 
dependent  environment  with  x  linearly  approximated  between  range  points:  (1)  the 
splitting  originally  proposed  by  Tappert  (Eq.  42),  has  a  local  error  (Eq.  55)  of 
order  A2  Independent  of  frequency;  (2)  the  later  splitting  (Eq.  44)  has  a  local 
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error  (Eq.  57)  of  order  A^  and  linear  in  frequency;  and  (3)  standard  Crank- 
Nicholson  (Eq.  60)  has  a  local  error  (Eq.  61)  of  order  A^  and  cubic  in  frequency. 
See  References  7,  14,  16,  17,  and  18  for  more  in-depth  analyses. 

Range-step-size  predictors  are  obviously  important  for  practical  implementation 
of  any  technique,  and  we  refer  the  reader  to  References  14  and  19  for  discussion  of 
some  of  the  techniques  that  are  utilized. 

The  advantage  of  split-step  over  Crank-Nicholson  is  that  it  is  more  efficient 
at  higher  frequencies  (linear  vs.  cubic)  permitting  a  larger  range  step  size.  Its 
disadvantage  is  that  a  uniform  z-grid  is  required,  and  in  situations  involving  rapi^ 
variations  in  depth  of,  for  example,  sound  speed,  a  large  transform  size  is  required 
(small  z-mesh  interval).  The  finite  difference  technique,  on  the  other  hand,  does 
not  require  a  uniform  grid,  allowing  for  the  flexibility  of  a  judicious  choice  of 
mesh  points  and  consequent  reduction  in  mesh  points. 

E.  HIGH  ANGLE  PE 


P 


I 


« 


1 .  General 


The  rational  linear  approximation  of  Section  G  is  not  appropriate  for 
solution  by  the  split-step  approach  because  of  the  variable  coefficients.  The 
techniques  employed  in  the  programs  presented  at  the  workshop  all  used  standard 
integration  in  r,  but  employed  special  representations  In  z  as  follows:  Greene  used 
cubic  splines  to  represent  the  z-dependence  for  his  rational  linear  approximation, 
whereas  Thomson  used  central  differences  and  Gilbert  used  finite  elements  for  the 
discrete  z-representation  of  Claerbout's  rational  linear  approximation.  Error 
analysis  in  range  Is,  in  general,  the  same  as  for  SPE  except  that  the  x  operator 
must  be  appropriately  redefined.  While  we  do  present  the  various  discrete 
techniques  for  representing  the  z-dependent  part  of  the  rational  linear  approxi¬ 
mation  we  do  not  include  a  corresponding  error  analysis.  Such  an  analysis  is  beyond 
the  scope  of  this  report. 

2.  Cubic  Splines 


Integration  over  r  of  Greene’s  Rational  Linear  Approximation,  Eq.  35,  using 
the  trapezoidal  rule  on  the  right  hand  side  and  neglecting  the  r  dependence  of  the 
index  of  refraction  leads  to 


(a 


A6  n  n+1 
2  '  v 


(y  -  * 


n  n+1 


zz 


(a  +  M)V  +  (y  +  -^)n*n 
2  2  zz 


(62) 


where  superscripts  n  and  n+1  refer  to  r*r„  and  rn+i,  respectively. 

In  order  to  treat  the  z  dependence,  the  following  spline  interpolation 
function  can  be  introduced  for 


v  + 


£  +  1 


v  + 


4Mr  (vj  -  v) 
£  i 
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(63) 
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where 


6 

i 


-  z 


l 


* 


v  =  (z  -  z4)/6a  , 


v  =  1  -  v  , 


(64) 


and  where  Is  the  evaluation  of  ^  at  z=»z^ .  The  are  parameters  to  be  deter¬ 
mined.  In  this  formulation  V'  and  4,zz  are  continuous,  and  the  conditions  on  ^at 
the  boundaries  z-zj  and  z-zn  require  gi=gn*0.  With  continuity  conditions 
for  and  p~*  $ z,  Eq.  63  yields  a  set  of  linear  equations  in  ^  and  which  can 
be  written  in  matrix  notation  as 

Mg  =  N  f  .  (65) 


Here  M  and  N  are  tridiagonal  and  g  and  ?  are  column  vectors  consisting  of  the  g 
and  ^ .  *■ 

Wow,  use  of  the  Interpolation  function  of  Eq.  63  in  Eq.  62  leads  to  an 
equation  of  the  form 


,n_n+l 
AL  * 


ARv  +  kftg n 


(66) 


where  an*'  are  diagonal  matrices  and  b[,  and  bg  are  constants. 

Finally,  multiplying  Eq.  66  by  M  and  using  Sq.  65  to  eliminate  Mg  we  arrive  at  the 
Implicit  form 


n  —  n+l 
CL' 


£3 


(67) 


This  equation  can  be  solved  numerically  in  an  efficient  manner  since  c£  and 
Cjj  are  tridiagonal.  Wote  that  this  formulation  assumes  p  “p(z),  and  results  in 
being  continuous  when  in  fact  It  should  be  discontinuous  for  discontinuous  p 
or  * nd*'s  of  refncrion.  However.  this  condition  can  be  remedied  with  a  slight 
modification  to  the  above  procedure  and  still  retain  the  tridiagonality , 

3.  Heterogeneous  Approximation 

Thomson  solves  Claerbout's  rational  linear  approximation  in  the  form 
appropriate  for  depth-dependent  density.  This  form  can  be  derived  froo  the 
Helmholtz  equation  directly  or  by  first  deriving  CPE  for  variable  density  and  then 
introducing  the  appropriate  expansion  similar  to  that  used  in  Section  Il-'t.  Without 
reproducing  the  derivation  here,  we  simply  state  the  result  that  GPE  forP»  p(z)  1$ 
llentlcal  to  that  given  by  Eq  .  8  except  that  q  19  redefined  as  q  where 


-  I  + 


kp  ~£z 


5s 


(68) 
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Likewise,  the  Claerbout  approximation  to  OPE  for  this  same  case  is  Identical  to  Eq. 
32  except  that  again  q  is  replaced  by  q.  Introduction  of  the  envelope  function 
then  leals  to  Claerbout's  rational  linear  approximation  for  variable  density: 

(1  +  ^  q  <P  .  (69) 

4  r  L 


Integration  of  Eq.  69  over  an  interval  rn+^~rn,  assuming  constant  coefficients 
over  the  interval,  leads  to 


(1  + 


!)n(i 
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(70) 


Standard  central  differencing  in  z  is  next  employed  to  evaluate  all  but  the  p depen¬ 
dent  part  of  q.  To  evaluate  this  Thomson  uses  the  "heterogeneous  approximation"20 
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(71) 


As  in  the  spline  formulation,  the  result  is  a  tridiagonal  system  of  equations  that 
is  easily  solved  numerically.  This  formulation,  however,  does  model  density  discon¬ 
tinuities  although  the  discontinuities  la  n  are  not  explicitly  treated. 


4.  Finite  Elements 

Gilbert  solves  Claerbout's  rational  linear  approximation  for  variable  den¬ 
sity  P(  z)  using  the  method  of  finite  elements.  The  discretized  equations  in  r  and  z 
were  actually  derived  by  using  a  variational  technique  to  minimize  the  appropriate 
Ldgrangtan  density.  However,  for  simplicity,  we  will  present  a  different  hut  equiv¬ 
alent  derivation  here. 


We  start  with  Claerbout's  rational  linear  approximation  for  variable  den¬ 
sity  p(z)  (Eq.  70).  The  z-dependence  is  now  treated  by  the  following  interpolation 
function: 

viz)  •  Y,  Cs4j(z)  ,  (72) 


where 
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In  this  expression  th»*  d>j  are  cubic  polynomials  within  an4  zero  oulside  a  z-griJ 
Interval,  rod  the  and  v.,  are  V’  and  at  the  z-grid  points. 

Insertion  of  the  representation  for  V  given  by  Eq.  72  into  Eq.  70,  multi¬ 
plication  by  integration  -jvej  z  Is* ds  to  a  linear  system  of 

equations  that  can  bo  expressed  In  matrix  fora  as 


(Hn  -  ik0AXniCn+!  =  [Mn  +  ik0AXrf:n 


nu) 
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where  M  and  N  are  block  tridiagonal  with  blocks,  and  C  is  the  column  vector  of 
Cj’s.  This  system  of  equations  can  oe  solved  numerically  in  an  efficient  manner, 
but  at  a  slower  rate  (for  the  3arae  z-mesh)  than  the  tridiagonal  systems  of  Greene 
and  Thomson.  However,  discontioui ties  in  both  p  and  n  are  treated  explicitly  and 
can  be  input  as  arbitrary  functions  of  z.  This  method  need  not  use  cubic  polynom¬ 
ials  for  <#> j ,  and  can  be  easily  generalized.  In  f«ct,  for  linear  d>j  the  result 
is  essentially  equivalent  to  using  second  order  differences  for  <f>  and  leads  to  a 
tridtagonal  system. 

F.  INITIAL  CONDITIONS  AND  WAVENUMBER  FILTERS 

As  stated  before,  the  principal  advantage  of  PE  over  elliptic  equation 
approaches  is  that  PE  allows  numerical  solutions  which  march  out  in  range  from  an 
initial  field.  However,  because  of  the  sensitivity  of  the  solutions  to  these  init¬ 
ial  data  and  rhe  inherent  far-field  nature  of  the  PE  approximations,  careful  spec¬ 
ification  of  the  depth  dependence  of  the  field  at  starting  range  rQ  is  required. 

Both  initially  and  down  range  it  is  advantageous  (and  sometimes  necessary)  to 
limit  the  size  of  the  vertical  wavenumbers  or  angles  of  propagating  energy.  The 
choice  of  limit  depends  on:  the  particular  version  of  PE  (regions  of  validity  of 
small  angle  uPE,  high  angle  PE,  etc.),  the  boundary  conditions  (especially  attenua¬ 
tion  of  high  angle  energy  by  the  bottom),  computational  constraints  (depth-mesh 
spacing  and  range  step  decrease  as  angular  aperture  increases),  and  application 
(e.g.,  study  of  narrow  beam  vertical  arrays). 

For  a  range-independent  environment  the  exact  initial  field  can  be  obtained  by 
a  normal  mode  solution.  In  a  range-dependent  environment,  especially  for  sources 
near  bottom  slopes,  it  is  noc  practical  to  obtain  the  exact  initial  field,  and  only  a 
number  of  approximate  approaches  are  available.  One  is  to  obtain  the  field  from  a 
r  ly  theory  calculation,  but  this  technique  nas  its  disadvantages  and  was  not  used  by 
any  of  the  PE's  exercised  at  the  workshop.  The  techniques  that  were  used  fell  into 
two  categories:  (1)  normal  modes  and  (2)  point-source  simulations,  each  of  which 
can  be  "filtered"  to  limit  vertical  wavenumbers.  These  approaches,  as  exercised  by 
models  for  the  workshop  test  cases,  are  discussed  in  order  below. 

The  normal  mode  startup  Is  straightforward.  Under  the  assumption  that  the 
ocean  Is  stratified  near  the  source,  the  normal  mode  solution,  v(r0,z),  to  the 
elliptic  wave  equation  for  the  initial  envi ronment  Is  calculated,  and  related  to  the 
initial  PE  field  ^(r0,z)  by: 

vlr0,z)  -  r,,*2  v(r0,7,)e"ik°ro  .  (75) 


The  value  of  rQ  Is  arbitrary  hut  Is  usually  on  the  order  of  several  wavelengths. 

Such  a  starter  can  of  course  he  used  at  a »/  range,  and  has  the  advantage  of  being 
e;act.  In  a  range-independent  context.  One  can  also  predetermine  the  tmaxiraum  angle 
(vertical  wive number)  of  propagating  energy  at  the  source  by  truncating  Che  normal 
mode  sum  at  the  appropriate  eigenvalue.  However,  It  may  not  always  be  necessary  to 
use  the  Initial  environment  to  obtain  the  normal  mode  solution  as  Is  demonstrated  by 
Gilbert  in  Sec.  VII  for  Test  Case  3,  where  he  greatly  simplified  the  initial  envi¬ 
ronment  ar*l  truncated  the  normal  mode  sura. 

Instead  of  depending  on  a  normal-mode  calculation  to  initialize  PE,  certain 
Implementation!’  use  an  analytical  initial  field  at  r0  designed  to  closely  match 
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the  spherical-wave  solution  of  the  elliptic  equation  for  a  point  source  at  depth 
zg  in  a  homogeneous  environment 


P(r,z) 


eik0R 

R ' 


R  =  (r^ 


+  (z  -  Zg)2)1* 


(76) 


This  approach  is  based  on  assumptions  that  neither  boundary  nor  refraction  effects 
are  important  to  range  rQ  and  that  the  far-field  approximation  obtains, 

V0»l. 

Following  this  logic  and  at  the  same  time  seeking  to  reduce  the  energy  of 
large-angle  propagation,  Tappert3  has  proposed  an  initial  PE  field  with  envelope 
of  form 


iKr0,z)  =  S  exp(-(z  -  zs)2/W2) 


(77) 


where  S  and  W  are  to  be  determined.  When  normalized,  this  function  of  z  has  the 
form  of  a  Gaussian  density  function  with  mean  zg  and  variance  W2/2,  and  is  thus 
termed  a  "Gaussian  source.”  If  not  truncated,  its  wavenumber  (kz)  power  spectrum 
(Fourier-transform  modulus  squared)  is  also  Gaussian,  with  mean  0  and  variance 
W2/4.  When  limited  to  t*,e  transform  has  sidelobes  caused  by  the  con¬ 

volution  of  the  Gaussian  with  a  shape  of  form  (siny)/y  where  y  is  (zmax‘kz). 

Following  Brock^,  the  solution  of  the  PE  equation  in  a  homogeneous  melium 
with  the  starting  field  given  by  Eq.  77  yields 


<P(r,z) 


S(1  *  1k^r’5exp<- 


(2  ~  2g) 

W2(l+  i  2r 

k^W7 


2 

-) 

•) 


(78) 


By  matching  the  elliptic  (Eq.  76)  and  parabolic  (Eq.  78)  solutions  for  small  angle 
(r  >>|z-zs|),  one  obtains 


1(2-)* 

Wk0 


and 


w  *» 


According  to  Tappert3  some  of  the  advantages  of  the  Gaussian  form  are: 

(1)  It  is  a  "natural"  for  a  snail  angle  approximation  (|z-z3|«r);  (2)  In  the 
limit  of  !z|»W,  it  ts  smooth  and  has  minimal  sidelobes  in  both  z  and  k  space; 

(3)  The  phase  la  correct  for  a  point  source  in  a  homogeneous  medium,  and  the  depth 
dependence  is  correct  to  “second  order  in  angle"  (i.e.  to  second  order  in 
|z-zsl  /O  • 


The  principal  disadvantages  of  the  Gaussian  field  were  pointed  out  by  D.  Wood 
and  J.  Hanna  In  talks  at  the  workshop:  (1)  It  is  symmetric  about  the  source-depth 
ami  hence  not  appropriate  near  boundaries;  (2)  In  wavenumber  (angle)  space  it  "rolls 
off"  slowly— *t hereby  allowing  energy  corresponding  to  high  angles  to  propagate  and 
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at  the  same  time  reducing  the  energy  for  small  angles  (i.e.,  it  does  not  match  the 
field  of  a  wavenumber-limited  point  source). 

There  are,  of  course,  alternatives  to  the  Gaussian  initial  field  which  also 
simulate  a  wavenumber-limited  point  source.  Thomson's  model  (this  workshop), 
employs  a  weighted  (sinx)/x  form  with  x  =  /?(z-z8),  whose  transform  (in  the  limit 
of  z»  /?  )  is  a  rectangular  "boxcar"  function  with  value  zero  for  |k|>/?  and  constant 
for  | k  j  <  /3 .  Garon,  Hanna,  and  Rost^l  and  Hanna  (this  workshop)  described  the 
source  function  for  the  SAI  PE  in  the  context  of  a  digital  (wavenumber)  filter.  The 
shape  of  the  inital  angular  (wavenumber)  spectrum  is  specified  in  a  digitized  form. 
The  transform  gives  the  initial  depth-dependent  field.  This  allows  for  special 
source  shapes,  and  can  easily  accommodate  changes  in  the  z-mesh  size.  Both  the 
Thomson  and  Hanna  approaches  give  a  sharper  cut  off  of  angles  than  the  Gaussian 
field — thereby  allowing  large  z-mesh  size  and  range  step. 

For  reasons  mentioned  above,  it  is  usually  desirable  to  limit  propagating 
angles  (wave  numbers)  not  only  at  the  source,  but  also  down  range.  Some  versions  of 
PE  (e.g.,  SAI)  accomplish  this  by  applying  a  wavenumber  filter  to  the  Fourier- 
transformed  solution  at  each  range  step.  This  is  especially  efficient  for  the 
split-step  algorithm  since  the  transform  is  available  at  each  step.  Without  such  a 
filter,  it  is  possible  for  PE  solutions  to  accumulate  contributions  from  what  are 
equivalent  to  modes  with  angles  outside  PE’s  realm  of  applicability. 

IV.  TEST  CASES  -  GENERAL 

A.  INTRODUCTION 

To  provide  a  basis  for  comparing  and  analyzing  PE  programs,  four  test  problems 
were  devised,  each  designed  to  test  some  aspect  of  the  program's  capabilities.  The 
tests  are  obviously  not  exhaustive  but  do  focus  on  some  of  the  capabilities  required 
to  meet  the  needs  of  the  SEAS  program.  The  choice  of  test  cases  was  subject  to  the 
following  constraints:  (1)  the  number  of  cases  should  be  small;  (2)  input 
parameters  that  define  the  problems  should  be  relatively  simple  (thereby  ruling  out 
many  range  dependent  environment  cases);  and  (3)  computer  costs  for  execution  should 
be  kept  as  low  as  possible  (thereby  ruling  out  long  range  and/or  high  frequency). 

The  four  problems  selected  were: 

Test  Case  1  -  Range  Dependent  Surface  Duct 

Test  Case  2  -  Bilinear  Profile 

Test  Case  3  -  Range-Independent  Shallow  Water  Environment 

Test  Case  4  -  Basin/Slope/Shelf  Configuration  for  Geoacoust leal  Bottom 

Details  of  these  test  cases  will  follow  in  later  sections. 

Participants  were  encourage!  to  exercise  their  PE  programs  for  as  many  test 
cases  as  possible  even  though  the  programs  may  not  have  been  designed  to  treat  ail 
cs  :s,  e-g..  ^plying  a  PE  program  to  Test  Case  3  that  Is  not  designed  to  treat 
attenuation,  density,  and  sound  speed  discontinuities.  The  objective  was  not  to 
encourage  a  poor  showing  but  rather  to  determine  the  consequences  of  applying  a  PE 
to  a  test  case  it  was  not  designed  to  handle  and  at  tne  same  time  determine  the 
consequences  of  applying  some  other  PE  that  was  designed  to  hand  le  that  same  t*»st 
case — thereby  permitting  an  assessment  of  the  value  such  modifications  or 
improvements.  That  objective  was  met  with  a  fair  degree  of  success. 
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Table  2  summarizes  tiie  characteristics  of  the  PE  programs  that  were  exercised 
against  the  test  cases.  In  this  table,  and  in  ones  that  will  follow,  we  use  the 
shorthand  notation  given  in  Table  3.  Abstracts  of  these  PE  programs  are  given  in 
Appendix  8.  Appendix  B  also  contains  abstracts  for  SNAP  and  FFP,  the  SACLANT  ASW 
Research  Centre  programs  that  were  used  as  "reference”  calculations.  The  SNAP 
result  is  used  as  the  reference  calculation  in  the  test  case  comparisons  and  is 
referred  to  in  the  figures  as  "MODE."  Table  4  shows  which  test  cases  were  attempted 
by  each  program.  The  X's  that  are  underlined  indicate  that  the  test  case  was  run 
but  the  results  are  not  included  in  this  report.  The  SAI-1  case  2A  was  not  to  the 
correct  scale,  there  was  a  yet  unidentified  problem  with  the  NOSC  results  for  case 
4,  and  the  SAI-2  program  contained  a  compiler  error  that  affected  all  but  cases  2B 
and  2C.  The  NUSC-2  program  is  not  a  typical  PE  (it  uses  a  normal  mode  program  to 
solve  SPE),  and  the  result  for  Case  2  is  included  primarily  for  comparison  with  the 
usual  PE  results. 

B.  COMMENTS  ON  WORKSHOP 

In  fairness  to  the  workshop  participants  some  general  comments  should  be  made. 
First,  there  were  scaling  problems  at  the  workshop  that  precluded  definitive 
comparison  of  results.  This  was  due  in  part  to  requiring  a  less-than-ideal  plot 
scale  for  viewgraph  presentation,  and  in  part  by  not  recognizing  the  consequences  of 
the  distortion  introduced  by  reproducing  viewgraphs  from  originals  on  a  variety  of 
machines.  Second,  some  participants  had  the  advantage  of  resources  (e.g.,  normal 
mode  programs)  with  which  they  could  "tune  up"  their  PE  results.  Ideally  two  types 
of  results  from  each  participant  would  have  been  extremely  useful:  (l)  a  routine 
run;  and  (2)  a  thorough  analysis  of  the  test  case  with  a  "tuned  up"  run.  Since 
resources  varied,  initial  results  spanned  this  range. 

Because  of  the  above  considerations,  it  was  decided  at  the  workshop  that  every¬ 
one  would  be  given  the  reference  calculations  (to  be  shown  in  later  sections),  and 
that  participants  could  then  submit  new  results,  due  a  month  after  the  workshop  (May 
1981),  based  on  this  information  and  knowledge  gained  at  the  workshop.  Participants 
were  also  offered  the  opportunity  to  comment  on  their  results,  in  writing,  to  aid  in 
the  analysis.  However,  the  response  from  the  authors  of  the  various  PE  programs  was 
minimal  and  the  interpretation  of  the  results  is  almost  entirely  due  to  the  authors 
of  this  report.  It  was  also  decided  to  include  the  results  of  Gilbert  (NORDA),  al- 
tluugh  lie  had  not  presented  them  at  the  workshop  because  the  program  was  described 
at  the  workshop,  the  results  were  delivered  before  May  1981,  and  the  results  pro¬ 
vided  additional  evidence  of  the  accuracy  of  the  "high-angle"  PE. 

One  additional  comment  is  appropriate  at  this  point.  Most  programs  exercised 
at  the  workshop  had  been  in  existence  for  some  time  and  were  likely  debugged.  Other 
programs  such  as  DREP  (Thomson),  NORDA  (Gilbert),  SAI-2  (Greene),  and  URI  (Papadak- 
ls)  were  relatively  new,  and  likely  were  not  thoroughly  debugged. 

C.  FORMAT  OF  TEST  CASE  RESULTS 

Standardization  and  conventions  are  obviously  necessary  for  comparison  of  test 
case  results.  To  this  end,  a  source  is  understood  to  be  the  fixed  field  point  at 
zero  range,  and  a  receiver  Is  Che  field  [joint  that  varies  with  range.  Unless  other¬ 
wise  noted,  the  sound  speed  is  assumed  to  extend  as  constant  beyond  the  deepest 
depth  specified,  i.e.,  an  Iso-veloclty  half  space  underlies  each  case  and  precludes 
the  return  of  energy  (a  radiation  boundary  condition  at  z-*-).  Tie  sound  speed  is 
assumed  to  be  a  linear  function  of  depth  between  given  [joints.  Transmission  loss  is 
in  dB  re  1  o  and  range  is  In  km.  Average  TL  results  are  intensity  averaged  over  a 
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Table  3.  Notation  used  in  test  cases 


SPE 

CIT 

HAC 

HAG 

G 

S 

D 

NM(r) 

U 

SS 

FD 

FE 

CS 

FULL 

IMP 

AB 

FB 

A 


? 


Standard  PE 

Tappert's  Co-independent  PE 
Claerbout's  high-angle  PE 
Greene's  high-angle  PE 
Gaussian  source 
SAI  source 
DREP  source 

Normal  mode  source  at  range  r 
User  input 

Split-Step  algorithm 

Finite  difference 

Finite  element 

Cubic  splines 

Full  bottom  treatment 

Impedance  boundary  conditions 

Attenuating  bottom,  smeared  sound  speed 

False  bottom  with  attenuating 
layer 

Range  increment 
Depth  Increment 
Not  relevant 
Not  known  to  authors 


Table  4.  Programs  exercised  against  test  ca  ;s 


TEST  CASES 


1 

2 

3 

4 

ABC 

ABC 

A  B 

A  B  C  D 

BTL 

XXX 

XXX 

X  X 

X  X 

DREP 

XXX 

XXX 

X  X 

X  X  X  X 

NORDA 

X  X 

NOSC 

XXX 

X 

NRL-1 

XXX 

XXX 

X  X 

X  X  X  X 

NRL-2 

XXX 

XXX 

X  X 

NRL-3 

X  X 

X  X  X  X 

NUSC-l 

XXX 

XXX 

X  X 

X  X  X  X 

NUSC-2 

XXX 

X  X 

SACLANT 

XXX 

XXX 

X  X 

X  X  X  X 

SAI-1 

XXX 

XXX 

X  X 

X  X  X  X 

SAI-2 

XXX 

X  X 

X  X  X  X 

URT 

XXX 

XXX 

X  X 

X  X  X  X 

REF 

XXX 

X  X 

XX  X 
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specified  interval.  It  is  worth  noting  that  participants  used  different  averaging 
schemes  and  it  should  be  recognized  that  this  can  contribute  to  differences  between 
results. 

D.  INTERPRETATION  OF  RESULTS 

In  the  following  sections  (V-VIII),  for  each  test  case:  (1)  the  test  problem 
is  defined,  (2)  the  PE  program  results  are  presented,  and  (3)  the  results  are 
Interpreted.  The  first  two  are  straightforward.  An  explanation  of  why  each  PE 
performed  as  it  did  is  not  so  straightforward,  and  is  best  left  to  the  author  of 
that  PR  since  he  alone  knows  how  he  parameterized  the  particular  test  case  to  fit 
within  the  limitations  of  his  PE.  It  would  be  a  monumental  undertaking  for  the 
authors  to  "dissect"  every  PE  for  every  test  case,  especially  Case  4.  The  analysis 
for  each  test  case  is  therefore  superficial  in  the  sense  that  we  speculate  on  the 
origin  of  differences  among  various  PE  results  based  on  treatments  of:  (l)  initial 
field,  (2)  interfaces,  (3)  radiation  conditions,  and  (4)  type  of  PE  equation.  Our 
Interpretive  remarks  tend  to  be  general  rather  than  specific,  and  are  by  no  means 
conclusive.  Finally,  evaluation  of  how  well  a  PE  did  is  somewhat  subjective  and 
depends  to  a  great  extent  on  what  features  are  the  most  important  for  the  case  in 
question.  The  results  speak  for  themselves,  and  the  reader  can  draw  his  own 
conclusions. 

V.  TEST  CASE  1  (Range-Dependent  Surface  Duct) 

A.  DEFINITION  OF  PROBLEM 

This  test  case  has  a  range-dependent  surface  duct  and  consists  of  profile  la 
from  0  to  20  km  in  the  first  region,  a  transition  from  la  to  lb  over  the  range 
Interval  20  to  30  km  In  the  second  region,  and  profile  lb  from  30  to  50  km  in  the 
third  region.  Pertinent  parameters  are  shorn  in  Table  5. 

The  profile  in  the  transition  interval  is  range  dependent  between  depths  of  0 
and  300  m  and  is  given  by: 

h(r)  -  300  -  10(r-20) , 

Ca(r)  -  1480  +  1.238(r-20), 

and 

Cd(r)  »  1500  +  0.571(r-20), 


where 


h(r)  is  the  depth  (m)  of  the  duct, 

Cs(r)  is  the  sound  speed  (m/s)  at  zero  depth, 
r  is  the  horizontal  range  In  km, 


and 


Cd(r)  is  the  sound  speed  (m/s)  at  the  bottom  of  the  duct. 

There  are  three  parts  to  this  test  case  as  defined  by  the  receiver  depths; 

(A)  25  m,  (8)  250  m,  and  (C)  400  ra. 

The  25  m  source  Is  In  the  duct,  a9  Is  the  25  m  receiver.  The  250  m  receiver  is 
initially  in  the  duct  in  the  first  region,  goes  from  in-duct  to  below-duct  in  the 
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Table  5.  Parameters  defining  Test  Case  1 


TEST  CASE  1 


Frequency  ■  25  Hz 
Source  Depth  *  25  m 


Receiver  Depths  »  25,  250  and  400 

m 

Max  Range  ■  50  km 

Profile  la: 

Depth  (m) 

Sound  Speed  (m/s) 

Density  (g/cm3) 

Attenuation  (dB/km) 

0 

1480 

1. 

0. 

300 

1500 

1. 

0. 

1000 

1460 

1. 

0. 

Profile  lb: 

Depth  (a) 

Sound  S^eed  (m/s) 

Density  (g/cm3) 

Attenuation  (dB/km) 

0 

1492.38 

1. 

0. 

200 

1505.71 

1. 

0. 

300 

1500.00 

1. 

0. 

1000 

1460 

1. 

0. 

Plots — TL  vertical  scale  60-130  dB  at  lOdB/lnch 
Separate  plots  for  each  receiver  depth 
0  to  50  km  at  5  km/ in 
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transition  region,  and  is  below-duct  in  the  last  region.  The  400  ra  receiver  i9 
always  below  the  duct. 

Modal  analysis  of  this  case  shows  only  one  trapped  (but  slightly  leaky)  mode  in 
the  first  region,  that  transitions  into  an  untrapped  mode  in  the  last  region.  This 
should  be  an  easy  test  for  any  PE  program  since  phase  errors  should  not  be  a  problem 
(there  is  only  one  mode  and  C0  can  be  selected  with  confidence),  since  the  range 
dependence  of  the  environment  is  weak,  and  since  attenuation,  sound  speed  and 
density  discontinuities  or  rapid  variations  are  not  present.  That  is  to  say  that 
all  PE  programs  should  agree  and  give  the  "correct"  answer.  See  Table  6  for  a 
summary  of  the  PE  program  parameters  used  for  Test  Case  1. 

B.  RESULTS  AND  SUMMARY  FOR  TEST  CASE  1 

Mo  reference  calculation  was  available  for  this  problem  over  the  full  range  (0 
to  50  km);  however,  the  test  results  of  NRL-1,  NUSC-1,  BTL,  and  SACLANT  are 
virtually  identical  and  agree  with  range- independent  environment  normal  mode  results 
up  to  20  km.  Figure  1  shows  the  NUSC-1  results  for  the  three  receiver  depths  and  is 
typical  of  NRL-1,  BTL,  and  SACLANT. 

The  results  from  other  programs  exhibited  slightly  different  behaviors,  and 
several  types  of  problems  were  uncovered  by  this  test  case.  These  problems  are 
examine!  using  part  C  (400  m  receiver)  as  the  example.  The  first  is  the  problem  of 
"false"  bottom  returns  and  Is  shown  in  Figure  2  comparing  NRL-2  with  NUSC-1  for  the 
400  m  receiver.  The  spurious  oscillations  in  the  NRL-2  result  are  due  to  medium-  or 
high-angle  energy  returning  into  the  problem  from  the  bottom,  i.e.,  caused  by 
insufficient  attenuation  in  the  bottom  or  by  transform  artifacts  (sidelobes  or 
aliasing)  which  redistribute  the  energy  into  higher  angles.  This  is  also  typical  of 
the  SAI-1  result  (not  shown).  Both  use  the  SAI  filter  (Hanna)  to  determine  the 
inltiaL  field  and  eliminate  higher  angle  energy  as  the  solution  progresses  in  range. 
The  difficulty  is  caused  by  the  algorithm  that  determines  the  mesh  size  (5)  and  prob¬ 
lem  extent  (zmax)  from  the  requested  beamwidth.  The  algorithm  yields  an  effective 
beamwidth  greater  than  or  equal  to  that  requested.  For  this  particular  problem  the 
hal f-beamwidth  is  greater  than  90°  regardless  of  the  beamwidth  requested,  ami  there 
is  little  or  no  attenuation  of  bottom  returns.  Brock  was  able  to  correct  the  NRL-2 
result  for  this  problem  by  using  a  new  zmax  inversely  proportional  to  the  beam- 
width.  The  result  (not  shown)  is  virtually  the  same  as  NUSC-1. 

False  bottom  returns  would  have  also  been  present  in  the  NUSC-1,  NRL-1,  BTL  and 
SACLANT  results  had  a  z;nax  =  W  been  used,  i.e.,  a  false  bottom  layer  of  333  m. 

The  problem  was  avoided  however  by  using  a  larger  layer  depth,  as  indicated  in  Table 
6.  The  degree  to  which  the  false  bottom  duplicates  the  radiation  condition  is 
obviously  dependent  on  the  attenuation  in  the  layer  and  the  layer  thickness. 

Oscillations  are  also  present  in  the  DRKP  result  compared  in  Figure  3  for  the 
400  m  receiver  depth.  Their  origin  is  not  easy  to  determine,  since  the  false  bottom 
used  in  th»»  DRKP  calculation  extended  to  6000  o.  The  spurious  oscillations  in  the 
URI  results,  shown  in  Figure  4  for  a  receiver  depth  of  400  m,  cannot  be  attribute*! 
to  a  false  bottom  since  impedance  conditions  were  employed.  Papadakis  attributes 
this  to  numerical  noise  that  can  occur  when  the  signal  gets  low;  however,  in  this 
case,  we  see  the  level  is  on  tV  order  of  90  dB  in  the  10  to  30  ko  interval. 

Finally,  the  NQSC  result  Is  shown  in  Figure  5  for  a  receiver  depth  of  400  ra. 

In  this  case  a  normal  mode  solution  was  used  to  generate  the  field  to  a  range  of 
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Figure  5.  NOSC  PE  results  for  Test  Case  1C  compared  with  NUSC-I  PE 


20  km  where  the  PE  calculation  is  initiated .  We  have  no  explanation  of  the  oscil¬ 
lations  between  20  and  25  km,  but  they  are  likely  due  to  matching  the  PE  calculation 
to  the  normal  mode  field. 

This  test  case  has  demonstrated  that  PE  results  for  a  narrow  beam  problem  with 
no  drastic  range  dependence  of  the  environment  and  no  bottom  complications  (other 
than  that  the  radiation  condition  be  satisfied)  are  relatively  Insensitive  to  com¬ 
puter  type,  starter,  z  and  r  mesh,  phase  corrections,  choice  of  CQ,  and  whether 
solved  by  SS  or  FD. 

Questions  raised  by  this  test  case,  but  not  answered,  are:  (1)  what  is  the 
most  effective  way  of  approximating  the  radiation  condition,  and  (2)  how  often 
should  a  range  dependent  environment  be  sampled? 

VI.  TEST  CASE  2  (Bilinear  Profile) 

A.  DEFINITION  OF  PROBLEM 

In  this  test  case  the  environment  is  range-independent  with  a  bilinear  sound 
speed  in  the  water  column  and  a  half-space  bottom.  It  consists  of  three  parts  and 
the  parameters  are  given  in  Table  7. 

The  rationale  for  this  case  is  to  test  the  beamwidth  capabilities  of  a  PE  pro¬ 
gram  without  the  complications  of  bottom  interaction  other  than  ensuring  that  no 
energy  returns  to  the  water  column.  Case  2A  requires  a  15°  half-beamwidth  capabili¬ 
ty,  2B  a  30°  half-beamwid* h  capability,  and  2C  a  40°  half-beamwidth  capability. 

There  is  a  continuous  spectrum  in  this  problem  because  of  the  branch  point  at  z=1500 
m;  however,  this  contribution  is  negligible  beyond  10  km  for  all  three  cases,  and 
SNAP  (which  does  not  treat  the  continuous  spectrum)  and  FFP  (which  does  treat  the 
continuous  spectrum)  should  agree,  as  they  did.  The  reference  calculation  for  this 
case  is  the  SNAP  result  which  we  refer  to  as  MODE. 

We  would  anticipate  that  all  PE  programs  should  be  able  to  do  well  on  2A  and, 
without  a  high-angle  capability,  they  should  have  some  difficulty  with  2B  and  have 
serious  problems  with  2C.  Initial  field  conditions  (amplitudes  cf  corresponding  PE 
modes)  and  phase  errors  will  govern  the  accuracy  of  the  results.  See  Table  8  for  a 
summary  of  PE  program  parameters  used  for  Test  Case  2. 

B.  RESULTS  FOR  CASE  2A  v 

For  this  test  case,  which  required  a  15°  half-beamwidth  capability,  all  PE's 
did  well  with  little  disagreement  in  detail.  The  most  accurate  result  was  NRL-1 
with  CPA.  We  do  not  show  it  because  it  is  virtually  the  same  as  the  MODE  result. 

We  show  instead  in  Figure  6  the  MODE  reference  calculation  along  with  the  envelope 
formed  by  the  DREP,  NRL-1  (no  CPA),  NUSC-1,  NUSC-2,  SACLANT,  and  URI  results.  This 
is  an  extremely  tight  envelope  and  there  is  little  to  choose  between  the  various 
results. 

Two  model  results  are  singled  out,  not  because  they  are  "bad",  but  because  they 
are  not  quite  as  accurate  as  the  other  and  are  to  some  degree  different.  The  first 
is  shown  in  Figure  7  which  compares  NRL-2  with  MODE.  The  "kinks"  in  the  NRL-2  curve 
are  due  to  straight-line  connections  between  range  point  solutions.  We  believe  the 
NRL-2  result  is  an  example  of  the  cumulative  errors  caused  by  a  crude  range  step 
size.  The  BTL  result,  shown  in  Figure  8,  is  not  so  easily  explained.  The  level  is 
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Table  7.  Parameters  defining  Test  Case  2 
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Figure  8.  BTL  PE  result  compared  with  reference  result  (MODE)  for  Test  Case  2A 
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good  but  we  are  unable  to  account  for  the  disagreement  in  detail.  We  would  have 
expected  this  result  to  fall  within  the  envelope  shown  In  Figure  6. 

One  conclusion  evident  from  this  test  case  is  that  CPA  is  effective  when  the 
SPE  phase  error  is  small,  seemingly  more  effective  than  CMOD  (used  in  the  SACLANT 
result  in  Fig.  6).  However,  SPE  without  corrections  is  very  accurate,  and  any 
corrections  should  be  small.  Since  this  is  the  only  case  where  CPA  is  used,  we 
cannot  evaluate  its  effectiveness  in  later  cases  where  the  SPE  phase  error  is  not 
small. 


C.  RESULTS  FOR  CASE  2B 

All  PE's  did  well  on  the  average  for  this  case  but  differed  in  detail.  A  half- 
beamwidth  capability  of  30°  is  required.  We  show  first  the  SAI-2  result  in  Figure  9, 
compared  with  the  MODE  solution.  This  is  the  most  accurate  PE  solution  for  this 
case  and  gives  surprising  agreement  in  detail  (the  pattern  in  Fig.  9  is  caused  by 
the  interference  of  22  modes).  This  is  an  excellent  demonstration  of  the  capability 
of  the  rational  linear  (high  angle)  approximation.  We  should  note  that  there  was  a 
compiler  error  in  Greene's  program  (discovered  after  the  workshop)  that  produced 
spurious  results  for  most  of  the  test  cases.  While  the  consequences  of  that  error 
are  included  in  Figure  9  (as  well  as  in  Fig.  13  for  Test  Case  2C) ,  we  assume  that 
elimination  of  the  error  can  only  improve  the  accuracy.  This  was  demonstrated  after 
the  workshop  by  results  Greene  obtained  using  a  some  what  different  program,  but 
still  a  rational  linear  approximation. 

In  this  same  connection  Figure  10  compares  the  DREP  result  against  the  MODE 
result.  Although  the  level  is  reasonable,  the  detailed  agreement  is  poor.  This  is 
difficult  to  account  for  since  the  DREP  program  also  solves  PE  in  a  rational  linear 
approximation. 

Results  from  the  remaining  programs  are  best  sorted  out  in  terras  of  whether  or 
not  the  CMOD  phase  correction  was  applied.  Figure  11  shows  the  envelope  formed  by 
the  NRL-1 ,  NUSC-1,  NUSC-2,  and  URI  results,  none  of  which  used  the  CMOD  correction 
while  Figure  12  shows  the  envelope  formed  by  the  BTL,  NRL-2,  SACLANT,  and  SAI-1 
results,  all  of  which  used  the  CMOD  correction.  These  figures  show  that  the  CMOD 
correction  produces  significant  differences,  greater  than  the  differences  within  the 
two  envelopes  caused  by  beamwidth,  mesh  size,  etc.  It  would  appear  that  CMOD 
improves  the  results. 

D.  RESULTS  FOR  CASE  2C 

Results  for  2C,  which  required  a  40°  half-beamwidth  capability,  are  similar  to 
those  for  Test  Case  2B.  Transmission  loss  levels  over  the  range  interval  are  fair, 
but  the  details  are  not  very  good  except  for  SAI-2  shown  in  Figure  13.  This  is 
strong  evidence  that  a  "high-angle"  PE  can  give  accurate  results  for  problems  having 
beamwldths  up  to  40°.  The  DREP  result,  Figure  14,  has  the  correct  level  but  again 
serious  phase  problems. 

As  in  Test  Case  2B,  the  remaining  results  are  best  presented  in  terms  of 
whether  or  not  a  CMOD  correction  was  applied.  Figure  15  shows  the  envelope  formed 
by  the  NRL-1,  NUSC-1,  NUSC-2,  and  URI  results,  none  of  which  used  the  CMOD  correc¬ 
tion,  while  Figure  16  shows  the  envelope  formed  by  the  BTL,  NRL-2,  SACLANT  and  SAI-1 
results,  all  of  which  used  the  CMOD  correction.  As  for  Tept  Case  2B,  the  differen¬ 
ces  between  the  envelopes  are  greater  than  the  differences  within  the  envelopes.  It 
is  not  clear  whether  CMOD  improves  the  overall  results  in  this  case. 
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- MODE 

- SAI-2 


TRANSMISSION  LOSS  (dB  re  1m)  _  TRANSMISSION  LOSS  <dB  re  1m) 


Figure  11.  Envelope  of  results  of  PE  programs  that  did  not  use  CMOD  compared 
with  reference  result  (MODE)  for  Test  Case  2B 


Figure  12.  Envelope  of  results  of  PE  programs  that  did  use  CMOD  compared  with 
reference  result  (MODE)  for  Test  Case  2B 
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TRANSMISSION  LOSS  (dB  re  1m) 


Figure  15.  Envelope  of  results  of  PE  programs  that  did  not  use  CMOD  compared 
with  reference  result  (MODE)  for  Test  Case  2C 


Figure  16.  Envelope  of  results  of  PE  programs  that  did  use  CMOD  compared  with 
reference  result  (MODE)  for  Test  Case  2C 
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E. 


SUMMARY  FOR  CASE  2 


Test  Case  2  has  demonstrated  that  SPE  can  give  an  accurate  answer  for  a  small 
beamwidth  problem  if  there  are  no  sound  speed,  density,  or  attenuation  disconti¬ 
nuities,  whereas  for  larger  beamwidth  problems  average  levels  are  predicted  well  but 
phase  errors  cause  disagreement  in  detail.  It  has  also  demonstrated  that  a  "high- 
angle"  PE,  in  this  case  Greene's  rational  linear  approximation,  can  provide  a  good 
approximation  to  the  elliptic  solution  for  a  range-independent  problem  containing 
half  beamwidths  up  to  40°.  In  addition  this  test  case  indicates  that  in  a  situation 
where  the  phase  errors  are  small  (i.e.,  SPE  is  very  accurate),  CPA  provides  a  more 
accurate  correction  than  CMOD.  For  larger  beamwidth  problems  all  SPE's  give  essen¬ 
tially  the  correct  level  but  are  unable  to  reproduce  the  details.  CMOD  appears  to 
slightly  improve  this  situation.  It  is  not  known  whether  CPA  would  provide  an 
improvement  in  general.  The  PE  results  are  relatively  independent  of  r-z  mesh, 
initial  field,  details  of  imposing  the  radiation  condition,  and  method  of  solution, 
i.e.,  split-step  or  finite  difference. 

VII.  TEST  CASE  3  (Shallow  Water) 

A.  DEFINITION  OF  PROBLEM 

The  environment  for  this  case  is  range-independent  and  consists  of  an  isoveloc¬ 
ity  water  column  over  an  isovelocity  half-space  bottom.  Parameters  for  the  two  parts 
of  the  problem  are  given  in  Table  9.  It  represents  an  idealized  shallow  water  envi¬ 
ronment  and  was  designed  to  test  whether  PE  can  accurately  handle  sound  speed,  den¬ 
sity,  and  attenuation  discontinuities  for  a  simple  problem  where  the  half-beamwidth 
requirements  are  less  than  20°. 

There  is  a  continuous  spectrum  in  this  problem  because  of  the  branch  point  at 
100  m.  However,  its  contribution  is  negligible  beyond  5  km,  and  thus  the  SNAP  and 
FFP  results  should  agree,  as  they  do  (comparison  not  shown).  There  are  11  propa¬ 
gating  modes  (corresponding  to  a  half-beamwidth  of  18.5°).  However,  with  the  source 
and  receiver  in  the  middle  of  the  water  column,  a  7  mode  solution  (corresponding  to 
a  12°  half-beamwidth)  is  fairly  accurate  in  the  5  to  10  km  range  and  a  9  mode  solu¬ 
tion  (corresponding  to  15°  half-beamwidth)  is  extremely  accurate  in  the  5  to  10  km 
range.  Thus  a  half-beamwidth  capability  between  12  and  15°  should  be  adequate  to 
accurately  handle  Case  3A.  For  Case  3B,  with  the  source  and  receiver  just  off  the 
bottom,  the  higher  modes  are  more  strongly  excited,  and  all  11  modes  are  required 
for  an  accurate  solution  in  the  5  to  10  km  range  interval.  Thi3  in  turn  requires  a 
half-beamwidth  capability  of  18.5°. 

It  should  be  noted  that  for  Case  3B,  where  the  source  is  within  0.08  wave 
lengths  of  the  bottom  interface,  we  also  have  a  test  of  the  initial  field  employed 
if  it  assumes  that  the  source  is  a  few  wavelengths  from  any  boundary.  See  Table  10 
for  a  summary  of  PE  program  parameters  used  for  Test  Case  3. 

As  mentioned  in  Sec.  III-F,  the  N0RDA  program  employed  as  a  starter  a  normal 
mode  solution  based  on  a  simplified  initial  environment.  This  environment  was  iso¬ 
velocity,  constant  density,  with  zero  absorption  and  with  pressure  release  boundary 
conditions  at  c  ■  0  and  z  »  650  m. 

B.  RESULTS  FOR  CASE  3A 

Results  for  this  case  are  mixed  and  for  the  first  time  some  PE's  have  incorrect 
levels.  It  is  worth  noting  that  CMOD  was  not  applied  for  this  te3t  case  (both  parts 
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Table  9.  Parameters  defining  Test  Case  3 

TEST  CASE  3 

Frequency  =  250  Hz 
Water  Depth  ■  100  m 
Sound  Speed  In  Water  3  1500  m/s 
Density  in  Water  *  1.0  g/cm^ 

Density  in  Bottom  »  1.2  g/cnr* 
Attenuation  in  Water  =  0 
Bottom  Attenuation  -  .5  db/A 
Bottom  Sound  Speed  *  1590  m/s 
Max  Range  ®  10  km 

Part  A: 

Source  Depth  -  50  m 
Receiver  Depth  »  50  m 

Part  B: 

Source  Depth  *  99.5  m 
Receiver  Depth  **  99.5  m 


Plots — Separate  plots  for  each  part  from 
5  to  10  km  at  0.5  km/ in. 


A  &  B)  as  it  is  not  appropriate  for  isovelocity  problems.  The  most  accurate  result 
was  that  of  NORDA  shown  in  Figure  17;  it  is  essentially  the  same  as  the  MODE  result. 
This  accuracy  is  the  result  of  proper  treatment  of  discontinuities  and,  as  will 
shortly  become  evident,  of  the  efficacy  of  the  rational  linear  approximation.  On 
the  other  hand  the  DREP  result,  shown  in  Figure  18,  seems  to  be  off  by  a  constant 
range  shift,  for  which  we  have  no  explanation.  The  compiler  error  encountered  by 
Greene  (SAI-2)  seriously  affected  the  resists  for  this  test  case  and  plots  are  not 
shown.  His  corrected  version  and  new  program  give  results  comparable  tc  the  NORDA 
results. 

The  next  test  results  are  shown  in  Figure  19,  which  displays  che  envelope  of 
the  NRL-1,  NRL-3,  NUSC-1,  SACLANT  and  URI  transmission  losses.  The  level  and  detail 
are  good,  and  the  envelope  is  fairly  tight  despite  the  difference  in  equations  being 
solved,  treatment  of  the  bottom,  starter,  and  solution  technique.  Notice  in  partic¬ 
ular  that  NRL-1  does  very  well,  although  it  does  not  have  a  bottom  treatment  appro¬ 
priate  for  this  problem  (see  Table  2). 

An  extremely  important  observation  can  be  made  upon  comparison  of  the  NORDA 
prediction  (Figure  17)  and  the  envelope  of  Figure  19:  the  high  angle  PE  is  more 
accurate  than  SPE  for  a  case  requiring  a  half-beamwidth  of  less  than  15°.  The  point 
made  by  this  result  is  that  high-angle  PE  is  not  only  more  accurate  in  phase  than 
SPE  for  large  angles,  but  is  also  more  accurate  in  phase  at  small  angles.  The  SPE 
phase  errors  in  this  test  case  are  evidently  just  enough  to  prevent  forming  the 
precise  mode  in;erference  pattern.  Since  phase  errors  are  accumulative  in  range, 
this  has  important  implications  for  relatively  narrow  beam,  long  range  propagation 
problems. 

The  NRL-2  result  is  presented  in  Figure  20.  This  result  is  not  bad,  but  is 
different  from  those  in  Figure  19.  The  difference, we  believe, is  primarily  due  to 
the  range  step  size  and  to  some  degree  the  treatment  of  the  bottom.  The  SAI-1  and 
BTL  results  are  shown  in  Figure  21.  Here  the  levels  are  considerably  in  error  and 
the  details  are  completely  different. 

C.  RESULTS  FOR  CASE  3L 

This  was  an  extremely  difficult  case  as  illustrated  by  the  comparison  of  MODE 
(SNAP)  and  FFP  in  Figure  22.  Agreement  is  excellent  except  near  the  null  in  the 
vicinity  of  7  ka.  There  are  obvious  numerical  problems  that  can  arise  when  11  modes 
must  almost  exactly  cancel.  Other  normal  mode  results  (not  shown  here)  Indicate 
that  the  FFP  solution  is  correct;  however,  we  will,  as  for  the  previous  lest  cases, 
use  the  MODE  (SNAP)  results  as  reference.  The  null  at  7  ka  Is  on  the  order  of  1/2 
km  wide  anti  is  Interesting  In  Its  own  right.  It  is  also  predicted  by  a  fully 
coherent  ray  theory  calculation.  Ray  theory  snaly  is  indicates  that  there  is  a 
complicated  Interaction  of  well-known  effects  which,  by  happenstance ,  combine  to 
produce  a  mill.  It  Is  not  appropriate  to  present  that  analysis  here,  nor  to  examine 
the  circumstances  under  which  It  can  occur.  Of  more  present  interest  is  that  the 
analysts  shows  that  Its  reproduction  requires  an  accurate  treatment  of  the  ;>hase  and 
water/bottom  interface. 

This  test  case  posed  problems  for  PE  both  In  1 *vel  and  in  detail;  however,  the 
NORDA  result  shown  In  Figure  23  is  very  accurate  (recall  that  the  MODE  result  d«>es 
not  faithfully  reproduce  the  null).  This  Is  further  supporting  evidence  that  a 
high-angle  PE  with  proper  bottom  treatment  fs  an  extremely  powerful  approximation. 
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Figure  17.  NORDA  PE  result  compared  with  reference  result  (MODE)  for 
Test  Case  3A 
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Figure  13.  DREP  PE  result  compared  with  reference  result  (MODE)  for  Test  Case  3A 
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Figure  24  shows  the  envelope  of  the  NRL-1,  NUSC-1,  and  URI  results,  all  SPE’s. 
Notice  that  the  average  level  is  good  but  the  detailed  pattern  is  not  good.  Figure 
25  is  the  SACLANT  result  and  is  shown  separately  as  an  indication  of  how  C0  can  be 
selected  to  reproduce  a  particular  feature.  In  this  case  the  null  at  7  km  is  repro¬ 
duced,  but  the  details  at  other  ranges  are  no  more  accurate  than  the  Figure  24 
results.  The  NRL-2  transmission  loss  is  shown  in  Figure  26;  it  is  neither  better 
nor  worse  than  the  SPE’s  in  Figure  24,  but  is  different.  The  DREP  and  NRL-3  results 
are  shown  in  Figure  27.  These  again  are  no  more  or  no  less  accurate  than  the 
results  in  Figures  24  and  26.  The  SAI-1  and  BTL  results  are  shown  in  Figure  28  and 
quantitatively  the  same  as  for  Test  Case  3A,  i.e.,  the  levels  and  details  are  off. 

Finally,  note  that  a  display  of  the  sensitivity  of  transmission  loss  to  source 
treatment  is  obscured  by  the  lack  of  agreement  of  any  method,  except  NORDA’s,  with 
the  reference.  Even  though  the  source  is  .08  wavelengths  from  the  water/bottom 
interface,  the  fact  that  the  NORDA  normal  mode  starter  is  based  on  a  homogeneous 
medium  (no  discontinuities)  suggests  that  the  form  of  the  initial  field  is  not  very 
important  for  this  case.  Other  results,  such  as  the  envelope  of  Figure  24,  in  which 
models  utilizing  a  variety  of  source  types  yield  similar  predictions,  support  this 
conclusion. 

D.  SUMMARY  FOR  CASE  3 

We  repeat  here  several  specific  results  of  this  test  case  which  appear  to  be 
generallzable.  First,  PE  can  do  well  in  predicting  average  transmission  loss  levels 
for  bottom  interaction  problems  only  if  some  reasonable  approach  to  handling  the 
bottom  is  included  in  the  PE  algorithm.  This  approach  may  be  crude  in  some  speci¬ 
fics,  such  as  in  this  case  where  several  programs  used  constant  density  even  though 
there  was  a  density  change.  The  conclusion  holds  whether  SPE  or  other  PE  forms  are 
used  and  whether  the  PE  is  implemented  through  split-step,  finite  difference,  or 
finite  element  techniques.  Second,  correct  detail  in  the  range  dependence  of  the 
transmission  loss  is  lacking  unless  energy  at  even  relatively  low  angles  (10°)  is 
treated  accurately.  Third,  the  high  angle  approximation  can  represent  this  energy 
much  more  accurately  than  SPE  and,  because  of  this,  can  track  transmission  loss 
through  a  larger  range  of  levels.  The  last  point,  one  which  has  long  been  accepted, 
is  that  the  choice  of  C0  appears  to  have  little  effect  on  average  level,  but  can 
have  major  effects  on  specific  detail,  and  in  fact  can  be  chosen  to  produce  "special 
effects"  such  as  the  proper  location  of  the  deep  fade  in  part  B. 

VIII.  TEST  CARE  4  (Basin/Slope/Shelf  with  Geoacoustical  Bottom) 

A.  DEFINITION  OF  PROBLEM 

This  test  case  has  a  range-independent  sound  speed  in  the  water  column,  a 
range-dependent  bottom  depth  consisting  of  a  basin-slope-shelf  configuration,  and  a 
simple  geoacoustic  bottom  (a  sediment  layer  overlying  a  basement)  with  sound  speed 
proportional  to  that  at  the  bottom  of  the  water  column.  The  parameters  defining 
this  case  are  given  in  Table  11.  This  is  perhaps  the  most  realistic  of  the  test 
cases  and  is  a  modified  version  of  a  recent  SEAS  problem  of  interest.  The  700  m 
receiver  has  a  positive  depth  excess  of  1346  m  in  the  basin,  becomes  bottom  limited 
on  the  slope  at  a  range  of  171  km,  and  strikes  the  bottom  slope  at  a  range  of  192 
km.  The  150  m  receiver  is  bottom  limited  throughout.  Because  of  the  test  case 
parameter  values,  energy  up  to  angles  of  40°  returns  to  the  water  column.  This  case 
is  therefore  a  stringent  test  for  any  PE  because  it  requires  proper  treatment  of 
interface  conditions  (discontinuous  sound  speed  and  density  at  the  water  column- 
sediment  and  sediment-basement  interfaces),  a  large  half-beamwidth  capability,  and 
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ENVELOPE:  NRL-1,  NUSC-1,  URI 
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Figure  24.  Envelope  of  NRL-1,  NUSC-1,  and  URI  PE  results  compared  with 
reference  result  (MODE)  for  Test  Case  3B 
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Figure  25.  SACLANT  PE  result  compared  with  reference  result  (MODE)  for 
Test  Case  3B 
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-  MODE 

- NRL-2 


Table  11.  Parameters  defining  Test  Case  4 


TEST  CASE  4 


Frequency  *  25  Hz,  Source  Depth  ■  600  m,  Maximum  Range  ”  250  km 


WATER  COLUMN  PARAMETERS 


Depth(m) 

C(m/s)  Density  (g/cm^) 

Attenuation  (dB/km) 

0 

1539.3  1. 

0. 

30 

1539.8  1. 

0. 

200 

1534.2  1. 

0. 

600 

1502.4  1. 

0. 

700 

1495.4  1. 

0. 

800 

1491.8  1. 

0. 

1000 

1488.0  1. 

0. 

1100 

1487.5  1. 

0. 

1200 

1487.9  1. 

0. 

3410 

1525.0  ‘  1. 

BOTTOM  DEPTH 

0. 

Range 

Bottom  Depth 

0-150  km 

Constant  3410  m 

150-200  km 

Linearly  decreasing  with  increasing  range  from  3410 
to  200  m 

200-250  km 

Constant  200  m 

SEDIMENT 

Sound  speed  at  top  of  sediment  is  CTa(r)  ■  .975 

*  Cy(r)  where  0y(r)  is 

the  sound  speed  at  the 

bottom  of  the  water  column. 

Sound  speed  at  the  bottom  < 

sediment  is  CBa(r)  ■  1. 

.305  x  CTs(r) 

Sediment  thickness  -  454  m  (constant  in  range) 

Attenuation  *  0.015  dB/km/Hz 
Density  »  1.5  grn/cm^ 

PLOTS:  700  a  RECEIVER  DEPTH 

4a  -  1.0  km  average  from  10  to  50  km  at  4  km/ in 
4b  -  20  km  average  from  0  to  250  km  at  25  km/ln 

150  a  RECEIVER  DEPTH 

4c  -  1.0  km  average  from  175  to  225  km  at  5  ka/in 
4d  -  20  km  average  from  0  to  250  km  at  25  ka/in 

For  additional  clarification  at  a  range  of  200  km:  Water  depth  *  200  a;  Cy  (200 
ka)  -  1534.2  m/a;  CB„  (200  km)  -  1952.1,  CTS  (200  ka)  -  1495.8;  and  Basement 
depth  *  654  a. 


a  range  sufficiently  great  (250  km)  that  cumulative  phase  errors  can  significantly 
distort  the  resultant  PE  field* 

Details  of  normal  mode  predictions  for  the  range-independent  portion  of  this 
problem  were  found  to  be  sensitive  to  particular  mode  program's  treatment  of  sound 
speed  interpolation  in  depth,  density  discontinuities,  etc*  It  was  therefore 
decided  to  avoid  comparisons  of  the  details  of  the  field.  It  was  decided  instead  to 
first  compare  intensity  average  results  on  a  fine  enough  scale  (1  km)  in  particular 
regions  of  interest  to  allow  meaningful  comparison  of  the  relative  behaviors,  yet 
not  so  fine  a  scale  as  to  make  comparisons  unduly  difficult,  if  not  meaningless.  At 
the  other  extreme,  a  gross  intensity  average  (20  km)  provides  an  overall  comparison 
of  levels  and  trends.  Since  this  20  km  average  is  roughly  equivalent  to  an 
incoherent  mode  sum,  the  PE  results  should  be  relatively  insensitive  to  phase 
errors. 


Another  and  perhaps  more  important  driving  factor  in  the  use  of  intensity 
averages  versus  detailed  transmission  loss  is  that  the  reference  calculation  (SNAP 
or  FFP)  can  be  used  with  confidence  only  in  the  basin  area  (0-150  km).  The 
"correct"  solution  in  the  slope  and  shelf  areas  (150-250  km)  is  not  known  and 
comparison  of  results  is  somewhat  of  a  "bootstrap*'  situation.  Although  an  adiabatic 
SNAP  solution  is  used  as  a  reference  result  for  the  150  m  source  (20  km  intensity 
average  case  from  0  to  250  km),  it  should  be  recognized  that  it  is  an  approximate 
solution  and  its  accuracy  relative  to  the  "true  solution"  is  not  presently  known. 

A  coupled  mode  solution  for  Test  Case  4  would  have  been  desirable  but  was  not 
available.  Thus,  we  will  be  comparing  PE  results  that  contain  full  coupling  but 
have  beamwidth  constraints  with  an  adiabatic  mode  result  which  has  no  beamwidth 
constraints  but  ignores  mode  coupling.  Consider  an  additional  complication  in 
comparing  average  results:  for  any  detailed  intensity-range  function,  the  average 
result  will  be  sensitive  to  the  weighting  of  the  values  within  the  given  interval 
(e.g.,  uniform  or  nonuniform)  and  to  the  range  values  for  which  the  average  is 
generated  (e.g.,  a  running  range  step  average,  or  an  average  for  every  range  value 
corresponding  to  the  mid-point  of  the  range  interval). 

The  test  rationale  for  each  of  the  cases  is  as  follows: 

4A.  Capability  to  treat,  on  a  1  km  average  level  of  detail,  a  case  that  has 
positive  depth  excess,  but  considerable  bottom  interaction. 

4B.  Capability  to  correctly  partition  water  column  and  bottom  interaction 

energy  for  a  case  that  has  positive  depth  excess  but  considerable  bottom 
interaction. 

4C.  Capability  to  track  on,  a  l  km  average  level  of  detail,  up  a  slope  and 
onto  a  shelf  under  strong  bottom  limiting  conditions. 

4d.  Capability  to  correctly  partition  energy  across  a  basin  up  a  slope  and 
onto  a  shelf  under  strong  bottom  limiting  conditions. 

For  Test  Ca9e  4  we  found  it  difficult  to  identify  patterns,  similarities, 
etc.,  among  the  various  PE  results  which  would  allow  sensible  groupings  for 
graphical  presentations.  Rather  than  presenting  the  results  individually  (wbi:h 
makes  intercooparison  almost  Impossible),  the  results  are  all  presented  in  one 
figure  for  each  test  case.  Table  12  gives  the  PE  program  parameters  for  Test  Case 
4.  For  convenience,  in  Table  12  we  use  the  following  short-hand  notation  to 
describe  the  type  of  Intensity  averaging  used: 
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LIN  -  equal  weighting  of  points  within  the  averaging  Interval 

EXP  -  exponential  weighting  of  points  within  the  averaging  interval 

R  -  average  over  the  interval  for  each  range-step 

I  -  average  over  the  interval  for  each  range  corresponding  to  the  mid-point 

of  the  interval . 

B.  RESULTS  FOR  CASE  4A 

The  PE  results  for  the  700  m  receiver  case  (1  km  intensity  average)  are  mixed 
and  are  shown  in  Figure  29  with  the  MODE  result  given  by  the  black  line  and  the 
various  PE's  color  coded  according  to  the  legend.  It  is  clear  from  Figure  29  that 
no  one  PE  does  well  throughout  the  range  from  10  to  50  km.  As  an  aid  in  understand¬ 
ing  these  results  see  Figure  30  which  shows  the  MODE  result  for  modes  1  through  16 
(corresponding  to  a  12.5°  half-beamwidth),  modes  1  through  56  (corresponding  to  a 

30°  half-beamwidth),  and  the  full  solution  modes  1  through  77  (corresponding  to  a 

40°  half-beamwidth).  The  16  mode  solution  contains  only  RR  paths  while  the  56  and 
77  MODE  solutions  include  bottom-interacting  paths.  Notice  the  small  spread  of  the 

various  PE  results  and  agreement  with  the  MODE  result  in  Figure  29  at  ranges  (34  to 

40  km  and  near  46  km)  where  in  Figure  30  the  RR  paths  are  important.  In  all  other 
regions  the  RR  paths  are  not  important,  and  the  spread  in  the  various  PE  results  is 
greater.  The  56  mode  solution  and  the  77  mode  solution  are  nearly  the  same  except 
in  the  range  interval  from  10  to  approximately  14  km.  Notice  also  in  Figure  29  that 
all  PE’s  over  estimate  the  loss  from  10  to  14  km  because  the  higher  angle  paths  are 
not  present.  Apart  from  the  RR  region,  detailed  agreement  is  poor  and  average 
agreement  is  fair>  Although  they  are  not  really  designed  to  handle  this  wide  a 
beamwidth,  the  NRL-1  and  SAI-1  results  compare  surprisingly  well;  this  is  probably 
due  to  the  large  half  beamwidths  (greater  than  30°).  The  primary  reason  for  lack  of 
better  agreement  of  PE  in  the  16  to  30  ka  region  is  the  impact  (on  a  1  km  intensity 
average  level  of  detail)  of  the  phase  errors  associated  with  SPE.  An  exception  to 

this  is  the  DREP  result  for  which  we  have  no  explanation. 

C.  RESULTS  FOR  CASE  4B 

A  reference  calculation  for  the  700  m  receiver  depth  (20  km  intensity 
average)  case  was  available  only  in  the  basin  region  and  Is  shown  in  Figure  31  as 
the  black  line  along  with  the  various  PE  program  results  (color-coded).  The  DREP 
result  appears  to  be  the  most  accurate.  Apart  from  the  NRL-1  result  which 
underestimated  the  loss  (and  which  employed  a  40°  half-beamwidth)  the  remaining  PE's 
tend  to  overestimate  the  loss  and  form  a  relatively  small  envelope  that  decreases  in 
width  with  increasing  range  to  less  than  2  dB  at  a  range  of  150  km.  This  envelope 
is  even  smaller  if  we  do  not  include  the  NUSC-1  and  URI  results  which  consistently 
overestimate  the  loss  on  the  order  of  a  dB  or  so  (probably  due  to  the  30° 
half-beamwidth  erployed).  In  the  plain  region  (0  to  150  km)  the  results  are  a 
measure  of  whether  a  PE  can  include  all  of  the  pertinent  energy  (up  to  40°)  and 
whether  the  interfaces  are  treated  accurately  enough  to  return  energy  from  the 
bottom. 


The  results  should  not  depend  strongly  on  phase  errors  since  the  20  km 
intensity  average  is  roughly  equivalent  to  an  incoherent  mode  sum.  Thus,  For 
programs  that  attempt  to  take  discontinuities  into  account,  the  transmission  loss 
should  be  greater  the  smaller  the  beamwidth.  This  appears  to  be  the  case. 

The  results  start  to  diverge  on  the  slope  (150  km)  to  the  range  where  the 
receiver  depth  and  bottom  depth  are  the  same  (192  km).  We  cannot  comment  on  the 
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Figure  29.  PE  program  results  (1  km  intensity  average)  compared  with  reference  result  (MODE)  for  basin 
of  Test  Case 


Figure  30.  Effect  of  numter  of  modes  {beam  width)  on  {1  km  intensity  average)  transmission  loss  in  Test 
Case  4A 


relative  accuracies  in  this  region  as  we  have  no  reference  calculation,  but  do  note 
that  the  spread  is  on  the  order  of  5  dB. 

D.  RESULTS  FOR  CASE  4C 

PE  results  for  Test  Case  4C  (150  m  receiver  depth)  are  shown  in  Figure  32, 
but  no  reference  calculation  is  available.  Information  obtained  after  the  workshop, 
but  not  included  here,  in  the  form  of  an  adiabatic  mode  calculation  (SNAP),  indi¬ 
cates  that  the  PAREQ  result  agrees  with  the  SNAP  result  in  general  level  but  not  in 
detail  (i.e.,  interference  pattern).  To  the  extent  of  our  confidence  in  SNAP,  we 
can  use  PAREQ  as  a  reference,  but  must  keep  in  mind  that  SPE  phase  errors  are  large 
enough  to  preclude  1  km  intensity  average  detail  comparisons. 

The  agreement  between  the  various  PE’s  in  Figure  32  is  generally  poor,  as 
perhaps  might  be  expected  because  of  the  different  treatments  of  interfaces,  method 
of  solution,  etc.,  but  the  large  degree  of  disagreement  is  a  little  surprising. 
Notice,  howe\ .  ,  that  there  tends  to  be  some  agreement  between  PE's  in  the  range 
Interval  182.5  to  187.5  km  (a  spread  of  2.5  dB  at  185  km).  This  apparently  is  a 
"slope-conversion"  effect,  caused  by  conversion  by  the  slope  of  small  angle  energy 
to  high  angle  energy.*  The  agreement  is  not  as  good  at  shorter  ranges  (7  to  10  dB), 
and  there  is  a  large  spread  (as  much  as  15  dB)  at  greater  ranges.  Also,  notice  that 
there  is  a  general  decay  in  level  beyond  200  km.  Although  the  level  is  different 
for  each  PE,  the  decay  rate  appears  to  be  the  same.  As  will  be  seen  in  Section  D, 
the  PE’s  enter  the  slope  region  (150  km)  with  a  spread  in  average  level  of  5  dB. 

This  spread  is  primarily  due,  just  as  for  the  600  m  receiver  case  (4A,  4B),  to 
initial  beamwidth  and  treatment  of  interfaces  and  radiation  boundary  conditions. 

Now  the  slope  is  approximately  3.7°  and,  therefore,  from  a  perspective 
each  encounter  with  the  slope  (water-column/ sediment  and/or  sediment/basement 
interface)  increases  the  angle  by  7.4°.  We  speculate  that  the  agreement  at  185  km 
is  due  to  (1)  I038  of  higher  angle  energy  that  contributed  to  the  disagreement  in 
the  plain  region  and  (2)  conversion  by  the  slope  of  the  RR  energy  (less  than  15°  in 
the  plain)  to  higher  angles.  However,  even  though  the  conversion  process  persists 
at  ranges  from  190  to  200  km,  the  various  PE's  are  now  each  converting  differently 
and  losing  the  higher  angle  energy  at  different  rates.  This  conversion  and  loss 
process  is  sensitive  to  the  treatment  of  the  slope  and  the  radiation  boundary 
conditions,  but  the  detailed  investigation  of  this  process  for  each  PE  is  beyond  the 
scope  of  this  report. 

E.  RESULTS  FOR  CASE  40 

A  reference  calculation  from  the  adiabatic  mode  program  SNAP  for  the  150  m 
receiver  depth  (20  km  intensity  average)  case  was  available  throughout  the  range  of 
the  problem  (0-250  km)  and  is  shown  in  Figure  33  a3  the  black  line  along  with  the 
various  PE  program  results  (color-coded).  It  is  interesting  to  note  that  the 
SACLANT  (PAREQ)  result,  which  was  accurate  on  a  1  km  intensity  average  basis  in  the 
slope-shelf  region  in  Case  4C  (175-225  km),  see  Figure  32,  displays  the  same 
accuracy  in  the  slope-shelf  region  on  a  20  km  intensity  average  basis  but  over¬ 
estimates  the  transmission  loss  in  the  basin  region  (0  to  150  km)  by  a  few  dli. 

*The  mechanics  of  this  process  are  complicated,  and  include  three  competing  effects 
on  transmission  loss:  (1)  increased  bottom  loss  as  paths  become  steeper,  (2)  de¬ 
creased  geometrical  spreading  loss  as  the  waveguide  constricts,  and  (3)  reduced 
surface  decoupling  (surface  image  interference)  loss  at  the  receiver  as  steep  paths 
become  more  important. 
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Figure  32.  PE  program  results  0  km  intensity  ayerage)  for  slope/shelf  of  Test  Case  4C 


Figure  33.  PE  program  results  ( 2^  knfr  intensity  average)  compared- to- reference  result  (MODE)  for  Test 
Case  4D 


Figure  33  is  a  dramatic  example  of  the  spread  in  average  level  (from  4  to  10  dB)  of 
PE  program  results  for  a  bottom  limited  wide  beamwidth  problem*  The  greatest  spread 
is  on  the  shelf,  where  we  have  the.  least  confidence  in  the  reference  calculation* 
Notice  that  in  the  basin  region  (0-150  km)  the  URI,  NUSC-1,  and  BTL  results  consist¬ 
ently  overestimate  the  transmission  loss  (all  use  a  half  beamwidth  of  308),  and 
both  URI  and  NUSC-1  have  a  reasonable  interface  treatment,  whereas  the  NRL-1  and 
SAI-1  results,  which  use  constant  density,  underestimate  the  transmission  loss 
(using  a  half  beamwidth  of  40°  and  35°,  respectively).  The  DREP  result  appears  to  be 
the  most  accurate  in  the  basin  and  SACLANT  on  the  slope  and  shelf  *  Again  as 
mentioned  in  part  C,  the  programs  with  smaller  beanwidth  (30°)  appear  to  have 
removed  too  much  energy  in  going  up  the  slope.  This  speculation  is  reinforced  if  one 
notes  that  the  decay  rate  is  about  the  same  for  most  of  the  PE’s  (and  MODE)  beyond 
about  205  km,  but  the  spread  in  level  is  larger  than  in  the  basin. 

F.  SUMMARY  FOR  TEST  CASE  4 

This  test  case  has  provided  some  interesting,  if  somewhat  disappointing,  re¬ 
sults.  In  general,  the  test  shows  the  importance  of  careful  treatment  of  interfaces 
in  problems  dominated  by  bottom-interaction  effects.  Also,  it  appears  that  a  high 
angle  PE  (appropriately  large  beamwidth)  is  needed  to  adequately  account  for  the 
bottom-interacting  energy,  especially  for  propagation  up  a  slope  and  onto  a  shelf. 
The  results  in  the  basin  region  appear  to  be  more  dependent  on  beamwidth  than  on 
type  of  PE  or  depth-mesh,  range-step  or  starter. 

Specifically,  these  PE  programs  perform  poorly  at  the  1  km  level  of  detail  in  a 
case  with  considerable  bottom  interaction  (Case  4A  and  4C).  In  tracking  up  a  slope 
and  onto  a  shelf  under  bottom  limiting  conditions,  they  are  not  consistent  (Case 
4C),  even  at  the  incoherent  (20  km  average)  level  of  detail  (Case  4D).  However, 
incoherent  levels  are  predicted  within  a  few  dB  in  the  range-independent  (basin) 
environment  (4B»  4D),  and  the  PE’s  do  very  well  even  on  the  small  scale  when  only  RR 
(low  angle)  energy  dominates  (4A  and  4C) . 

Test  Case  4  has  indicated  a  need  for  reference  solutions  to  a  set  of  range- 
dependent  environment  problems  that  could  serve  as  benchmark  test  cases. 

IX.  SUMMARY  AND  CONCLUSIONS 


This  report  has  presented  the  results  of  the  PE  workshop  held  from  31  March  to 
3  April  1981.  Two  major  objectives  were  met.  The  first  objective  was  to  provide  a 
forum  for  exchange  and  stimulation  of  ideas  for  those  presently  active  in  theoreti¬ 
cal  and  applied  PE  development,  as  related  to  underwater  acoustic  problems.  This 
objective  was  accomplished  partly  in  the  new  development  section  of  the  workshop  and 
partly  in  the  test  problem  section  of  the  workshop.  The  appendix  presents  the 
abstracts  for  the  new  development  portion  of  the  workshop  and  also  provides  brief 
model  descriptions  for  those  programs  exercised  for  the  workshop.  The  second  objec¬ 
tive,  a  comparison  of  PE  models  and  their  predictions,  was  covered  in  the  main  text 
and  is  discussed  below. 

The  first  half  of  the  main  text  reviewed  model  forms  and  solution  techniques 
for  PE  as  applied  to  underwater  acoustic  problems.  First  a  general  PE  and  then 
specific  PE  equations  were  discussed:  from  SPE  (historically  the  first  and  the  most 
common  form  in  use)  to  C0-independent  PE  and  rational  linear  PE  (the  latest 
development  in  general  use).  The  CMOD  and  CPA  modifications,  designed  to  correct 
SPE  phase  errors,  were  also  outlined.  Finally,  several  methods  of  solution  were 


discussed:  split-step  and  finite  difference  for  SPE,  and  schemes  for  the  three  high 
angle  PE’s  presented  at  the  workshop. 

The  second  half  of  the  main  text  presented  and  compared  PE  model  results  for  a 
set  of  four  ocean  acoustic  test  cases.  The  test  cases  ranged  from  a  simple  problem, 
expected  to  be  "easy"  for  any  PE  program,  to  a  fairly  realistic  basin-slope-shelf 
configuration  which  had  considerable  bottom  interaction  and  was  expected  to  be  a 
"difficult  test"  for  all  PE  programs . 

Detailed  comparisons  of  test  case  results  were  not  available  at  the  time  of  the 
workshop,  but  were  made  by  the  authors  afterward.  Conclusions  about  the  state  of  PE 
capability  were  then  drawn,  and  are  summarized  ?s  follows: 

1.  General 


a.  There  are  many  PE  models  in  operation  today.  They  differ  in  form,  in 
equation  solved,  in  treatment  of  inputs,  in  boundary  conditions  and  specifications, 
in  numerical  method,  in  computer  implementation,  and  in  display  of  the  output.  Each 
model  has  its  regime  of  applicability,  inherent  assumptions,  and  idiosyncrasies. 

Not  surprisingly,  this  diversity  of  models  was,  repeatedly,  exemplified  by  a 
diversity  of  test  case  results.  There  is  no  "standard,"  all-purpose  PE  model. 

b.  Nonetheless,  within  its  realm  of  applicability,  a  PE  model  can  be  a 
very  powerful,  efficient,  and  accurate  predictor  of  sound  transmission  properties. 
This  was  borne  out  in  those  test  cases  for  which  there  was  an  accurate  reference 
solution.  The  lack  of  a  reference  solution  for  certain  range-dependent  environment 
problems  serves  as  an  illustration  of  the  gap  in  technology  which  PE  can  fill. 

2.  SPE 


a.  The  test  problems  reconfirmed  that  the  accuracy  of  SPE  is  greatly 
dependent  on  the  vertical-angle  aperture  (beamwidth)  of  dominant  propagation  paths. 
From  a  modal  point  of  view,  the  smaller  the  extent  of  vertical  wavenumbers  or  phase 
speeds  for  propagating  modes,  the  smaller  the  SPE  errors  in  phase.  While  a  20° 
half-beamwidth  about  the  horizontal  (or  its  equivalent  "small  angle”  or  vertical 
wavenumber  limit  about  any  depression  angle)  is  a  traditional  rule  of  thumb  for 
tolerable  phase  errors,  in  SPE,  the  detailed  TL  prediction  suffers  from  the  effect 
of  cumulative  phase  errors  as  range  increases,  both  for  this  and  even  for  smaller 
apertures.  This  again  was  observed  in  the  test  problems. 

b.  For  problems  with  a  large  beamwidth  of  propagating  paths,  SPE  can 
predict  range-averaged  TL  (because  phase  errors  tend  to  average  out),  but  not  the 
details.  For  this  case,  and  to  a  lesser  extent  the  small  beamwidth  case,  the 
detailed  SPE  predictions  can  be  quite  sensitive  to  choice  of  C0. 

c.  The  accuracy  of  SPE  was  also  found  to  be  very  sensitive  to  the  treat¬ 
ment  of  interface  conditions,  especially  for  larger  beamwidth  problems.  Proper 
methods  of  dealing  with  abrupt  changes  of  sound  speed  or  density  with  depth  are 
required . 


d.  In  addition  to  the  inherent  limit  of  any  PE  to  "slowly"  range-varying 
environments,  SPE  was  found  to  be  sensitive  to  the  treatment  of  range  dependencies 
in  sound  speed  and  boundary  conditions.  This  is  as  much  an  enviroranental  as  an 
acoustic  problem. 
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e.  Test  case  results  tend  to  be  insensitive  to  choice  of  z  mesh,  range 
step  size,  solution  technique  (split-step  or  finite  difference),  and  starter, 
provided  of  course  that  all  were  within  reason. 

3.  CMOD 


The  CMOD  modification  to  SPT  was  found  to  reduce  SPE  phase  errors  and 
minimize  sensitivity  to  choice  of  C0,  within  the  limits  of  the  small  beamwidth 
constraint. 

4.  CPA 


The  CPA  modification  to  SPE  was  demonstrated  in  only  one  case,  where  SPE 
errors  were  small,  and  found  to  improve  the  results. 

5.  CIT 

The  C0-inde pendent  version  of  PE  showed  no  improvement  over  SPE  in  the 
test  problems. 

6.  High  Angle  PE 


Of  the  three  models  presented,  one  (NORDA)  did  only  Test  Case  3,  another 
(SAI)  had  results  only  for  Test  Case  2,  and  the  third  (DREP)  had  inconsistent 
results  (possibly  caused  by  implementation  problems  in  a  new  program).  General 
conclusions  cannot  be  drawn  from  this  small  sample,  but  individual  cases  suggest 
that  a  high  angle  PE  can  be  quite  accurate  for  beamwidths  as  large  as  +40°.  In 
addition,  the  finite  difference  and  finite  element  approaches  allow  interface 
conditions  to  be  accurately  represented. 

7.  General  Recommendations 


a.  There  Is  a  clear  need  for  bench  mark,  reference  solutions  to  a  set  of 
range-dependent  environment  problems. 

b.  PE  development  will  continue.  Persistent  evaluation  of  new  (and 
existing)  models  should  be  an  integral  part  of  that  development. 
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ALTERNATIVE  DATA  FOR  PE  MODELS 


David  H.  Wood 
Code  3342 

Naval  Underwater  Systems  Center 
New  London,  CT  06320 


The  Parabolic  equation  is  solved  by  numerical  models  that  go  forward 
in  range  one  step  at  a  time  using  previously  computed  values.  The  very 
first  values,  or  initial  data,  must  be  provided  to  the  model.  By  consider¬ 
ing  any  initial  data  whatsoever  after  expansion  as  a  weighted  sum  of  norm¬ 
al  modes,  we  can  immediately  see  what  the  solution  of  the  PE  would  be  if 
the  sound  speed  were  not  to  depend  on  range.  Matching  this  as  well  as  we 
can*  to  the  normal  mode  sum  that  satisfies  the  Helmholtz  equation,  and 
then  extrapolating  back  to  the  initial  range,  we  see  that  E0ntzo)0n(z)  is 
the  appropriate  initial  data.  There  are  three  deficiencies  in  the  result¬ 
ing  approximation  to  the  solution  of  the  Helmholtz  equation.  The  least 
important  Is  an  amplitude  error  for  each  mode  that  results  from  approxima¬ 
ting  the  square  root  of  the  horizontal  wave  number  for  each  mode  by  a 
single  representative  value*.  More  Important  is  that  each  mode  by  is  given 
a  phase  proportional  to  k(l-l/2sin29n)  instead  of  k  cos  9n  (the  hori¬ 
zontal  wave  number  of  the  n^  mode).  This  is  the  so-called  phase  error 
of  the  parabolic  approximation  and  is  small  for  small  angles  0n.  Most 
important  is  a  bizarre  difficulty  concerning  the  effervescent  modes  of  the 
Helmholtz  equation  that  decay  rapidly  with  range  because  cos  9n  is 
essentially  pure  imaginary.  For  ^hese  same  modes,  1-1/2  sin20n  is  essen¬ 
tially  pure  real  so  that  these  modes  propagate  in  the  parabolic  equation 
without  the  desired  exponential  decay. 


*Prof.  Tappert  pointed  out  during  the  presentation  that  the  initial  data 
0nUo°0n(z)  cos  “^On  would  give  an  improved  match  by  eliminating  ampli¬ 
tude  errors. 
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CONSTRUCTION  OF  A  BANDLIMITED  WAVENUMBER  SOURCE 
FOR  THE  PARABOLIC  EQUATION 

John  Hanna 

Science  Applications,  Inc. 

P.  0.  Box  1303 
McLean,  VA  22102 


All  solutions  to  the  Parabolic  Equation  (PE)  require  specification  of 
the  initial  pressure  field  as  a  function  of  depth  (or  wavenumber).  It  can 
be  shown  that  an  appropriate  specification  of  the  initial  field  in  verti¬ 
cal  wavenumber  space  for  beamwidths  of  30°  or  less  is  just 

0(o,k)  *  constant  for  -kb  £  k  <  kb 


«  0  otherwise, 


where 


kb  *  ko  sin0b, 
ko  -  w/c, 

9b  =>  desired  source  beamwidth. 

However,  since  most  numerical  solutions  to  the  PE  use  FFT's,  it  is  neces¬ 
sary  to  replace  this  rectangular  function  with  one  which  has  essentially 
the  same  shape,  but  also  acceptable  transform  properties.  For  the  SAI  PE 
model,  a  low-pass  finite- impulse  digital  filter  was  selected.  At  a  given 
frequency  this  filter  is  implemented  with  the  appropriate  number  of  points 
to  Insure  that  a  beamwidth  of  at  least  k0sin9b  is  achieved  la  wave- 
number  space.  The  principal  advantage  of  this  source  function  over  the 
traditional  Gaussian  function  is  that  it  minimizes  the  transform  size 
required  for  any  desired  beamwidth.  (The  difference  in  required  transform 
size  for  the  digital  filter  versus  the  Gaussian  is  typically  a  factor  of 
four.) 

If  one  constructs  the  ray  theory  equations  for  both  the  elliptic  and 
parabolic  equations  and  looks  for  t lie  equivalent  ray  trajectory  determined 
by  each  equation,  the  angle  of  a  ray  governed  by  the  elliptic  equation, 
9e,  is  related  to  Che  angle  under  the  parabolic  equation,  9p,  by 


sin  0e  •  fan  8p. 


So,  for  example,  the  90°  elliptic  equation  ray  is  equivalent  to  the  45° 
parabolic  equation  ray.  If  a  pressure  field  from  the  elliptic  equation  has 
wavenumbers  which  correspond  to  sin  9  >  1,  (i.e.,  k  >  k.G)  the  energy 
corresponding  to  them  is  exponentially  damped  out  as  the  field  propagates. 
However,  under  the  parabolic  equation,  energy  at  these  wavenumbers  is  wot 
damped  out  and  the  PE  equivalent  of  elliptic  rays  for  sin  8  >  1  will 
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propagate.  Such  energy  creates  unphysical  artifacts  in  the  PE  field.  To 
insure  against  these  artifacts,  the  digital  filter  wavenumber  cutoff  is 
constrained  to  the  following  interval 


k0sin  0^  £  k  <  k0 


which  achieves  the  desired  beamwiduh  and  prevents  the  artifacts. 

As  the  field  from  a  source  propagates  with  range  it  becomes  redis¬ 
tributed  in  depth  and  the  influence  of  refraction  and  bathymetry  is  to 
redistribute  the  energy  in  vertical  wavenumber  space.  Again,  because  most 
solutions  to  the  PE  are  based  on  discrete  Fourier  transfoims,  the  problem 
being  solved  is  implicitly  periodic  in  depth  and  wavenumber.  To  insure 
that  aliasing  does  not  occur  in  range  as  the  energy  redistributes  in  bot.ii 
depth  and  wavenumber,  the  low-pass  filter  is  also  used  at  each  range  step 
to  suppress  the  energy  at  high  wavenumbers  and  greatest  depths.  This  is 
accomplished  by  simply  multiplying  the  depth  and  wavenumber  fields  by  ^he 
filter  response  which  maintains  suitably  limited  depth  and  wavenumber 
fields . 

Early  portions  of  this  work  are  contained  in  a  talk  entitled  "Con¬ 
struction  of  a  New  Source  Function  for  the  Parabolic  Equation  Algorithm" 
by  H.  M.  Garon,  J.  S.  Hanna  and  P<  V.  Rost  (J.  Acoust •  Soc.  Am.  61  (51), 
512(A),  1977). 
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PULSE  RESPONSE  OF  PROPAGATION  CHANNELS  NEAR  CAUSTICS 


Lan  Nghiem-Phu 
Daubin  Systems  Corporation 
104  Crandon  Boulevard  #315 
Key  Blscayne,  FL  33149 

and 

Frederick  D.  Tappert 
University  of  Miami 
RSKAS 

4600  Rickenbacker  (SVY) 

Miami,  FL  33149 

A  Pulsed  Parabolic  Equation  (PPF)  model  has  been  developed  for  k',.e 
purpose  of  analyzing  the  time-domain  response  of  propagation  channels  for 
coherent  broadband  sources.  The  method  is  based  on  multi-frequency  PE  runs 
and  Fourier  transforms  of  the  complex  field  amplitudes  to  yield  the  ana¬ 
lytic  signal  in  the  time  domain.  The  time  domain  signals  are  normalized 
such  that  transmission  loss  (re  1  yd)  becomes  a  function  of  time  and  may 
be  resolved  into  separate  arrivals.  The  relative  phase  is  also  obtained  as 
a  function  of  time  and  is  nearly  stationary  through  the  center  of  single 
arrival . 

A  careful  study  of  the  accuracy  of  the  PPE  model  has  been  carried  out 
in  the  time  domain  by  comparison  to  ray  trace  calculation  of  arrival  times 
and  transmission  losses  of  single  arrivals.  In  the  neighborhood  of  caus¬ 
tics  for  pairs  of  arrivals,  we  have  obtained  some  new  results  on  the  pulse 
shapes. 

A  detailed  study  has  been  made  for  a  range  independent  environment 
consisting  of  a  bi-linear  profile  (iso-velocity  ocean  overlying  a  constant 
gradient  loss-less  bottom)  and  a  directional  source  having  200  Hz  band¬ 
width  and  200  Hz  carrier  (100  Hz  <  f  <  300  Hz).  The  pulse  response  was 
examined  at  several  ranges  as  a  function  of  depth.  Ir.  particular,  at  15  km 
range,  a  three-arrival  region  is  obtained  by  ray  theory  with  two  bounding 
caustics,  outside  of  which  there  are  single-arrival  regions. 

Pulse  response  amplitudes  were  compared  to  known  exact  solutions 
(spherical  spreading)  for  the  direct  arrivals  in  the  ieo-velocity  region. 
After  correcting  for  source  bo.ira  pattern  effects,  the  PPE  transmission 
loss  was  found  to  be  accurate  to  within  0.5  GB  out  to  15  ka. 

The  accuracy  oi  the  PPE  arrival  times  derived  from  phase  information 
was  found  to  depend  somewhat  on  the  choice  of  reference  sound  speed  Cp 
and  also  on  the  equivalent  ray  angle  with  respect  to  the  horizontal.  Ft  a 
poor  choice  of  CQ  and  large  angle  (about  24*),  the  error  in  travel  tlccs 
to  a  range  of  15  ka  could  be  as  much  as  10  msec.  However,  tor  an  optimal 
choice  of  C0  and  for  angles  less  than  IB*,  we  found  that  the  crrcr  in 
travel  times  ro  the  same  range  was  reduced  to  less  than  1  msec. 
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It  was  observed  that  the  peak  of  the  amplitude  was  not  situated  at 
the  geometrical  caustic  location  but  was  displaced  towards  the  illuminated 
side.  There  was  significant  energy  on  the  dark  side  which  decreased  expo¬ 
nentially  away  from  the  geometric  caustic  boundary,  and  resulted  in  pulse 
spreading  due  to  the  effective  low-pass  filter  caused  by  diffraction  in 
the  shadow  region.  Complicated  interference  patterns  developed  as  the 
caustic  was  approached  from  the  illuminated  side.  Separated  arrivals 
merged  together  when  their  relative  phases  differed  by  an  integral  mul¬ 
tiple  of  a  cycle  of  the  carrier  frequency.  For  arrivals  that  have  inter¬ 
acted  with  the  caustic,  there  was  evidence  of  pulse  deformation,  resulting 
In  a  trailing  tail. 

Similar  studies  are  in  progress  for  a  complex  environment  with  range 
dependent  profiles,  complicated  bathymetry,  and  range  and  depth  dependent 
horizontal  currents. 


THREE-DIMENSIONAL  ACOUSTIC  PROPAGATION 
US  INC  THE  PARABOLIC-EQUATION  METHOD 

J.  S.  Perkins,  R.  N.  Baer,  E.  B.  Wright,  and  R.  L.  Dlcus 

Code  5160 

Naval  Research  Laboratory 
Washington,  DC  20375 


We  have  extended  the  parabolic-equation  method  for  acoustic  propaga¬ 
tion  to  three-dimensions.  A  split-step  technique  results  in  an  algorithm 
for  marching  the  solution  in  range.  We  also  present  an  approximation  to 
this  full  three-dimensional  (3D)  algorithm  which  we  call  the  Nx2D  method. 
In  this  approach  we  solve  N  two-dimensional  problems  and  combine  the  re¬ 
sults  to  form  an  approximate  three-dimensional  solution.  This  technique 
allows  wave-front  bending,  which  is  the  dominant  mechanism  for  water-borne 
propagation. 

Wc  show  in  several  examples  how  horizontal  variations  in  the  ocean 
environment  affect  the  spatial  distribution  of  energy  and  the  distribution 
of  energy  In  azimuthal  angle  as  seen  by  an  array  of  hydrophones  placed  in 
the  calculated  complex-valued  acoustic  field.  In  the  first  example  we  con¬ 
sider  a  very  strong  linear  cross-range  gradient,  designed  to  produce  an 
obvious  and  predictable  effect.  Results  from  the  Nx2D  method  are  in  excel¬ 
lent  agreement  with  the  results  produced  using  the  3D  algorithm.  In  anoth¬ 
er  example,  we  investigate  propagation  through  a  strong  (but  realistic) 
ocean  eddy.  The  presence  of  the  eddy  between  source  and  receiver  causes  a 
change  in  the  apparent  source  bearing  of  0.5s.  In  this  case  the  Nx2D  meth¬ 
od  again  proves  to  bo  an  excellent  approximation.  In  the  final  example,  we 
use  the  full  3D  algorithm  to  simulate  the  performance  of  a  bottom-mounted 
array  of  hydrophones.  The  results  suggest  that  reflection  and  diffraction 
from  bathymetric  features  may  have  a  significant  effect  on  the  performance 
of  bottomed  arrays. 
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NUMERICAL  CALCULATIONS  OF  BOUNDS  ON  ARRAY  PERFORMANCE 
USING  THE  PARABOLIC-EQUATION  TECHNIQUE 

David  R.  Palmer 

National  Oceanic  and  Atmospheric  Administraton 
Atlantic  Oceanographic  and  Meteorological  Labs 
Ocean  Acoustics  Lab 
4301  Rickenbacker  Causeway 
Miami,  FL  33149 

and 

Marilyn  L.  Blodgett 
Code  5120 

Naval  Research  Laboratory 
Washington,  DC  20375 


In  order  to  evaluate  the  feasibility  of  proposed  under-sea  surveil¬ 
lance  systems  it  is  necessary  to  know  the  conditions  under  which  the 
effects  of  ocean  variability,  i.e.,  internal  waves,  on  system  performance 
can  be  ignored.  The  standard  approach  to  this  question  is  to  measure  or 
calculate  the  mutual  coherence  function  for  the  acoustic  field  since  it 
enters  directly  into  the  calculation  of  measures  of  system  performance 
such  as  array  gain.  Because  of  the  great  difficulties  in  determining,  both 
experimentally  and  theoretically,  the  mutual  coherence  function,  we  have 
taken  a  different  point  of  view.  Rather  than  attempting  to  calculate, 
e.g.,  array  gain,  we  seek  only  to  bound  it.  Consequently,  we  calculate  a 
limit  on  the  degradation  in  system  performance.  The  advantage  of  this 
approach  is  that  the  limit  depends  only  on  the  mean  pressure  field,  a 
quantity  far  easier  to  calculate  than  the  mutual  coherence  function.  The 
calculation  of  the  mean  field  does  require,  however,  the  use  of  a  general¬ 
ization  of  the  Markov  approximation  appropriate  when  diffracting  and 
refracting  effects  are  present.  With  this  approximation  the  calculation  of 
the  mean  field  can  be  carried  out  with  a  slightly  modified  version  of  the 
parabolic  equation  technique  of  Tappert  and  Hardin.  A  case  study  is  made 
for  propagation  in  the  Western  North  Atlantic  with  the  following  results: 
Independent  of  the  length  of  a  horizontal  array,  its  gain  will  not  be 
reduced  by  more  than  3  dB  for  propagation  to  l  convergence  zone  at  260  Hz, 
4  convergence  zones  at  130  Hz,  and  at  least  15  convergence  zones  at  50  Hz. 
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PRELIMINARY  VERSION  OF  A  FINITE-DIFFERENCE  PE  MODEL 


Ding  Lee  and  George  Botseas 
Code  3342 

Naval  Underwater  Systems  Center 
New  London,  CT  06320 

A  general  purpose  computer  model  which  calculates  propagation  loss, 
based  on  an  implicit  finite-difference  solution  to  the  parabolic  wave 
equation,  will  soon  be  available.  This  model  is  designed  to  offer  a  number 
of  useful  features:  modifications  are  easy  to  incorporate;  it  is  easy  to 
use;  and  reasonably  accurate  results  can  be  obtained  with  reasonable 
computational  speed.  In  shallow  or  deep  water  propagation,  the  present 
model  is  capable  of  handling  multiple  interfaces,  irregular  bottoms, 
arbitrary  bottom  boundary  conditions,  and  artificial  bottoms.  All  test 
results  show  excellent  agreement  with  published  results.  This  model  was 
used  to  solve  a  set  of  NORDA  workshop  problems.  These  results  will  be 
discussed . 


HIGH  ANGLE  PE 

Robert  R.  Greene 
Science  Applications,  Inc. 

1710  Goodridge  Dr. 

P.  0.  Box  1303 
McLeon,  VA  22102 

High  Angle  PE  is  a  full-wave  range  dependent  propagation  model 
similar  in  concent  to  PE.  It  is  designed  to  handle  the  dual  problem  of 
bottom-interac  to  ind  high  angles  of  propagation  up  to  about  45  degrees. 
It  is  based  upon  an  implicit  tridiagonal  range  step  and  a  cubic  spline 
discretization  in  the  vertical  coordinate. 

The  most  important  observation  about  the  method  is  that  since  it  uses 
a  finite  difference  approach  to  the  problem  rather  than  the  split-step  al¬ 
gorithm,  a  far  better  approximation  to  high  angle  propagation  is  used  in 
the  form  of  a  "rational  parabolic"  as  opposed  to  "parabolic"  approxima- 
t  ion . 


The  method  also  allows  an  accurate  treatment  of  strong  gradients  in 
sound  speed  and  discontinuities  in  sound  speed  and  density.  These  are  the 
main  difficulties  for  PE  in  modeling  bottom  interaction.  Strong  gradients 
are  not  a  limitation  for  a  finite  difference  technique  because  the  solu¬ 
tion  does  not  depend  on  splitting  an  operator  into  nearly  commuting  parts. 
Discontinuities  in  sound  speed  and  density  are  treated  accurately  by 
matching  cubic  splines  at  the  interfaces.  The  boundary  conditions  of  the 
wave  equation  are  the  natural  matching  conditions  for  the  splines. 

The  major  application  of  the  method  will  be  over  regions  where  bottom 
interaction  and  range  dependent  bathymetry  are  driving  the  problem.  Howev¬ 
er,  it  can  be  used  to  model  long-range  propagation.  It  should  be  particu¬ 
larly  appropriate  for  continental  slope  and  shallow  water  applications. 

The  work  is  currently  being  funded  by  0NR  Code  486. 
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A  FINITE  ELEMENT  METHOD  FOR  THE  PARABOLIC  WAV?  EQUATION 


Kenneth  E.  Gilbert 
Code  320 

Naval  Ocean  Research  and  Development  Activity 
NSTL  Station,  MS  39529 


A  semi-discrete  elliptic  wave  equation  is  derived  and  used  to  develop 
a  parabolic  equation  which  allows  propagation  angles  of  up  to  40°  with 
respect  to  the  horizontal.  The  derivation  is  based  on  the  minimization  of 
a  Lagrangian  in  which  the  vertical  dependence  of  the  acoustic  field  is 
writter  explicitly  in  terms  of  finite  elements  and  the  nodal  values  of  the 
pressure  3nd  vertical  particle  velocity.  Since  the  resulting  wave  equation 
is  in  terms  of  continuous  field  variables,  discontinuities  in  sound  speed 
and  density  are  easily  handled. 


A  FINITE-DIFFERENCE  TREATMENT  OF  INTERFACE  CONDITIONS 

Ding  Lee 
Code  3342 

Naval  Underwater  Systems  Center 
New  London,  CT  06320 

and 

Suzanne  T.  McDaniel 
Applied  Research  Laboratory 
The  Pennsylvania  State  University 
University  Park,  PA  16802 


A  finite-difference  approach  is  introduced  to  handle  both  horizontal 
and  Irregular  interface  conditions.  These  interface  conditions  are  devel¬ 
oped  to  preserve  continuity  of  pressure  as  well  as  continuity  of  vertical 
particle  velocity.  The  formulation  of  interfaces  is  by  means  of  a  finite- 
difference  scheme  in  conjunction  with  the  Taylor  series  expansion,  while 
the  handling  of  the  irregular  interfaces  is  by  line  segmentation.  The  com¬ 
plete  mathematical  development  will  be  outlined.  The  validity  of  this 
approach  is  clearly  shown  by  the  theory,  however,  test  results  will  be 
discussed  as  supporting  evidence. 
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A  HYBRID  PARABOLIC  EQUATION  CODE 


H.  Brock 
Code  1721.2 

Naval  Research  Laboratory 
Washington,  DC  20375 

The  Tappert-Hardin  split-step  Fourier  algorithm  is  an  accurate  and 
efficient  numerical  technique  for  the  solution  of  the  parabolic  wave 
equation,  particularly  on  a  vector  computer.  Difficulties  arise,  however, 
in  simulating  boundary  conditions  at  interfaces.  Discontinuities  introduce 
complications,  and  the  range  step  required  is  strongly  dependent  on  the 
gradient  of  the  index  of  refraction.  Such  boundary  conditions  can  be 
handled  by  finite  difference  techniques,  although,  in  general,  for  a  given 
accuracy  many  more  mesh  points  are  required  in  both  dimensions.  In  addi¬ 
tion,  on  a  vector  computer  the  techniques  available  for  the  solution  of 
the  finite  difference  equations  tend  to  be  less  efficient  than  the  split- 
step  algorithm.  Since  in  many  problems  RSR  propagation  dominates  over 
substantial  ranges,  on  a  vector  supercomputer  it  is  feasible  to  combine 
the  split-step  Fourier  algorithm  in  RSR  dominated  regions  with  a  finite 
difference  algorithm  in  regions  where  bottom  interactions  are  important. 
The  resulting  hybrid  PE  code  changes  both  the  parabolic  approximation  and 
the  solution  technique  used  dynamically  with  range  during  the  execution  of 
a  given  problem. 


EQUIVALENT  BOTTOM  FOR  USE  WITH  THE 
SPLIT-STEP  PARABOLIC  EQUATION  ALGORITHM 

Homer  Bucker 
Code  531 

Naval  Ocean  Systems  Center 
San  Diego,  CA  92152 

The  split-step  parabolic  equation  algorithms  (SSPEA)  have  become  the 
standard  Navy  model  for  acoustic  calculations  in  range  dependent  ocean 
areas.  This  occurred  because  SSPEA  is  extremely  simple,  runs  rather  quick¬ 
ly  on  a  computer,  and  does  a  good  job  of  modeling  the  refraction  of  sound 

caused  by  changes  of  sound  speed  with  depth  and  range.  However,  SSPEA  has 

problems  when  used  in  areas  of  strong  bottom  interaction.  These  are  due  to 
density  discontinuities  at  the  ocean  bottom,  sediment  rigidity  which  re¬ 
sults  in  shear  waves  in  bottom  sediments,  and  sound  scattering  from  ocean 
bottom  irregularities. 

The  equivalent  bottom  is  a  set  of  absorbing  liquid  layers  that  have 
the  same  density  as  the  bottom  water  and  thickness  equal  to  the  vertical 
sampling  interval  used  in  the  SSPEA.  The  sound  speed  and  density  of  each 

layer  is  adjusted  so  that  the  equivalent  bottom  has  the  “correct"  bottom 

reflection  properties.  That  is,  the  equivalent  bottom  should  have  the  same 
plane  wave  reflection  coefficient  as  a  function  of  grazing  angle  as  the 
actual  ocean  bottom.  In  addition,  restrictions  on  the  maximum  change  in 
sound  speed  and  attenuation  between  adjacent  layers  insures  that  range 
steps  in  the  SSPEA  can  be  kept  as  large  as  possible.  Also,  insofar  as 
scattering  effects  can  be  incorporated  into  the  reflection  coefficient, 
they  are  automatically  transferred  to  the  sound  level  calculation. 
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A  HIGH  FREQUENCY  RAY  BASED  PARABOLIC  EQUATION 


Frederick  D.  Tapper t 
University  of  Miami 
RSMAS 

4600  Rickenbacker  (SNY) • 
Miami,  FL  33149 


No  abstract  was  submitted  for  this  presentation. 


IMPEDANCE  FORMULATION  OF  THE  BOTTOM  BOUNDARY  CONDITION 
FOR  THE  PARABOLIC  EQUATION  MODEL  IN  UNDERWATER  ACOUSTICS 


John  .5.  Papadakls 
Department  of  Mathematics 
University  of  Rhode  Island 
Kingston,  RI  02881 


In  Implementing  Impedance  conditions  along  the  bottom  Interface  for 
the  parabolic  equation,  it  is  assumed  that  the  subbottom  region  is  homo¬ 
geneous,  semi-infinite  and  the  field  satisfies  a  radiation  condition  at 


The  derived  Impedance  is  a  condition  on  the  total  parabolic  field 
p(r,ag),  aQd  in  the  form  of  an  integral  equation  along  the  bottom  inter¬ 
face  Z-Zg 

r  ko  o 

U  /*  i-=-  (n  2  -  l )(r  -  s) 

I  pB  /  2  B 

p(r*ZB)  =  “  [ik7p’)^  ’  /  G  '(r  '  S>  PZ(s*zB)ds  ♦ 

o 

where  Pg/P  are  the  densities  in  the  bottom/water,  k»2lTf/c0  and  ng  is 
the  index  of  refraction  in  the  bottom. 

The  Integral  equation  is  approximated  by  a  finite  sum  and  is  incor¬ 
porated  in  an  implicit  finite  difference  algorithm  for  the  parabolic 
equation.  The  resulting  model  advances  the  parabolic  field  only  in  the 
water  column  and  does  not  require  the  introduction  of  an  artificial  sub¬ 
bottom  layer. 

Research  is  in  progress  on  models  with  sloping  bottom  interfaces  as 
well  as  elastic  subbottom  effects. 
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PERUSE:  A  PE  CODE  INCLUDING  ROUGH  SURFACE  SCATTERING 


Levis  B.  Dozier 
Science  Applications,  Inc. 
1710  Goodridge  Dr. 

P.  0.  Box  1303 
McLean,  VA  22102 


A  new  model  has  recently  been  developed  which  treats  rough  surface 
scattering  in  the  context  of  PE.  The  surface  is  assumed  to  be  piecewise 
linear  and  frozen  in  time.  A  sequence  of  conformal  mappings  locally  flat¬ 
tens  each  linear  surface  segment,  whereupon  the  split-step  Fourier  algo¬ 
rithm  can  advance  the  solution  for  one  step  in  each  transformed  space. 

Numerical  PERUSE  (PE  Rough  Surface)  results  for  surface  ducts  are 
analyzed  by  resolving  them  in  depth  into  Labianca's  virtual  modes  ("mode- 
matching’*).  The  resulting  modal  intensities  are  then  plotted  as  a  function 
of  range,  and  decay  rates  are  estimated  for  each  mode.  Comparison  with 
simple  Rayleigh  theory  shows  no  obvious  correlation.  A  fullwave  coupled 
mode  model  is  currently  being  investigated  to  explain  the  PERUSE  results. 


PROPAGATION  MODELING  WITH  THE  PARABOLIC  EQUATION 


F.  B.  Jensen  and  W.  A.  Kuperman 
SACLANT  ASW  Research  Centre 
La  Spezia,  Italy 


A  aeries  of  basic  propagation  problems  are  solved  using  the  parabolic 
equation.  First  we  show  results  of  some  mode  coupling  studies*-^  com¬ 
prising  propagation  in  a  wedge-shaped  ocean  (up-slope  and  down-slope)  and 
propagation  in  a  shallow-water  channel  with  an  abrupt  profile  change. 
Then  we  apply  the  PE  model  to  study  Lloyd-mirror  patterns.  These  results 
are  interesting  since  we  here  obtain  field  plots  that  explicitly  show  the 
effect  of  the  small-angle  approximation  inherent  in  the  PE  technique.  In 
this  connection  we  also  study  the  diffraction  of  Lloyd-mirror  beams  over 
the  top  of  a  wedge.  Finally,  we  demonstrate  the  application  of  PE  theory 
to  modeling  sound  propagation  over  a  seamount.  Here,  theory  is  compared 
with  broadband  experimental  data. 

1.  F.  B.  Jensen  and  W.  A.  Kuperman,  "Sound  propagation  in  a  wedge-shaped 
ocean  with  a  penetrable  bottom,"  J.  Acoust.  Soc.  Am.  67_,  1364-1566  (1980). 

2.  F.  B.  Jensen  and  W.  A.  Kuperman,  "Range-dependent  bottom-limited  propa¬ 
gation  modeling  with  the  parabolic  equation",  in  Bottom-Interacting  Ocean 
Acoustics,  edited  by  W.  A.  Kuperman  and  F.  B.  Jensen  (Plenum  Press,  New 
York,  1980)  pp.  451-466. 


84 


APPENDIX  B.  PE  PROGRAM  ABSTRACTS 


BTL: 

J.  D.  Seals 

86 

DREP: 

D.  J.  Thomson 

87 

NORDA: 

K.  E.  Gilbert 

87 

NOSC: 

D.  Gordon 

88 

NRL-1: 

J.  S.  Perkins  and  R.  N.  Baer 

89 

NRL-2/3: 

H.  Brock 

90 

NUSC-1: 

D.  Lee  and  G.  Botseas 

91 

NUSC-2: 

David  H.  Wood 

92 

SACLANT: 

F.  B.  Jensen  and  H.  R.  Krol 
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and  R.  Stieglitz 
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Robert  R.  Greene 
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John  S.  Papadakis 
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SNAP: 

F.  B.  Jensen  and  M.  C.  Ferla 
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FFP: 

H.  W.  Kutschale 

98 

BTL 

ARRAY  PROCESSOR  PE 


J.  D.  Seals 

Bell  Telephone  Laboratory 
Whippany  Road 
Whippany,  NJ  07981 


The  BTL  version  of  the  parabolic  code  (APPE)  was  developed  to  provide 
production  runs  with  turnaround  times  on  the  order  of  minutes  rather  than 
hours.  This  was  achieved  by  careful  implementation  of  the  code  on  an  AP- 
190L  array  processor  interfaced  to  an  UN1VAC  1108  host.  The  full  potential 
of  the  array  processor  was  realized  only  by  wholly  implementing  the  P.E. 
approximation  and  step  size  algorithm  in  the  array  processor.  While  the 
array  processor  propagates  the  code  forward  one  or  more  steps,  the  host 
sets  up  physical  input  for  the  next  step(s).  When  the  array  processor 
finishes  the  propagation,  a  minimal  I/O  exchange  occurs  and  the  process 
repeats  until  completion.  A  pre-cursor  routine  looks  ahead  to  warn  of 
major  changes  in  physical  properties  so  that  step  reduction  occurs  before 
rather  than  after  the  fact. 

Additional  features  of  the  BTL  code  include  both  linear  and  triangu¬ 
lar  sectors  interpolation  of  velocity  profiles,  range  dependent  Q-^  and 
bottom  sound  speed  gradients,  pseudo  discontinuities  at  the  water/sediment 
Interface,  mixed  data  units,  full  diagnostic  capabilities,  data  checking, 
linear  and  exponential  averaging  of  TL  data,  and  extensive  graphic 
capabilities. 

Although  the  current  version  of  the  code  was  developed  for  produc¬ 
tion,  its  highly  modular  and  structured  format  provides  an  excellent  base 
for  future  program  study  and  development. 
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DREP 


A  FINITE-DIFFERENCE  PE  MODEL 
D.  J.  Thomson 

Defence  Research  Establishment  Pacific 
FMO  Victoria,  B.  C. 

CANADA  VOS  1B0 

A  finite-difference  scheme  (hastily  invoked  for  the  NORDA  workshop) 
is  described  for  the  so-called  heterogeneous  formulation^  of  the  para¬ 
bolic  equation,  which  allows  for  spatial  variation  in  the  material  prop¬ 
erties.  The  sound-speed,  density  and  absorption  may  be  specified  at  each 
grid  point  of  a  finite-difference  mesh,  and  the  (internal)  boundary  con¬ 
ditions  are  everywhere  satisfied  implicitely.  Discontinuities  in  material 
properties  can  be  accommodated.  The  Crank-Nicolson  algorithm  is  used  to 
march  out  the  numerical  solution  to  the  parabolic  equation-  A  source 
function,  designed  to  be  band-limited  in  vertical  wavenumber,  is  used  to 
begin  the  solution  process.  Promising  features  cf  this  approach  include 
the  ability  to  treat  parabolic  equations  with  variable  coefficients,  the 
use  of  non-uniform  step  sizes  in  depth,  and  the  possibility  of  adapting 
recently  developed  absorbing  boundary  conditions^  to  limit  the  vertical 
extent  of  the  computational  grid.  While  not  thoroughly  tested  (or  even 
possibly  debugged),  the  new  code  was  applied  to  each  of  the  suggested  test 
problems. 

R.  Kelly,  R.  W.  Ward,  S.  Treitel,  and  R.  M.  Alford.  Geophysics,  Vol. 
41(1),  pp.  2-27,  (1976). 

^R.  W.  Clayton  and  B.  Engquist.  Geophysics,  Vol.  45(5),  pp.  895-904, 
(1980). 


NORDA 

A  FINITE  ELEMENT  HIGH  ANGLE  PE 

K.  E.  Gilbert 
Code  120 

Naval  Ocean  Research  and  Development  Activity 
NSTL  Station,  MS  39529 

The  NORM  program  solves  the  wave  equation  in  Cl.o,»rbout ' s  rational 
linear  approximation,  i  .e .  a  high-angle  parabolic  equal. on.  The  discrete 
representation  is  based  on  the  minimization  oF  a  Lagrangian  in  which  the 
vertical  dependence  cf  the  acoustic  field  is  written  explicitly  in  terras 
ol  finite  elements  and  the  nodul  values  of  the  pressure  and  vertical  par¬ 
ticle  velocity.  Since  the  resulting  wave  equation  is  in  terras  of  continu¬ 
ous  field  variables,  discontinuities  in  sound  speed,  attenuation  and 
density  are  easily  handled.  This  program  presently  treats  only  range  inde¬ 
pendent  environments.  A  normal  mode  solution  is  used  for  the  initial 
field,  and  the  radiation  condition  at  the  bottom  is  approximated  by  a 
"false"  bottom.  The  program  is  Implemented  on  the  NORDA  Cyber  171 
computer. 
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NOSC 


PARABOLIC  EQUATION  COMPUTER  PROGRAM 
D.  Gordon 

Naval  Ocean  Systems  Center 
San  Diego,  CA  92152 


This  underwater  acoustic  propagation  loss  program  s  a  split-step 
parabolic  equation  algorithm  using  an  FFT  with  no  corrections  or  exten¬ 
sions.  It  is  designed  to  be  started  from  a  normal  mode  field  input  at  any 
desired  range.  Sound  speed  profiles  are  normally  entered  as  linear  layers 
in  1/C2.  A  constant  attenuation  can  be  added  to  any  layer.  The  attenua- 
ating  layer  has  a  large  value  of  attenuation  at  the  bottom  of  the  attenu¬ 
ating  layer,  and  decays  exponentially  in  dB/Ar  upwards  so  that  attenuation 
is  negligible  at  the  top  of  the  attenuating  layer  and  above  it.  The  oper¬ 
ator  specifies  Az,  Ar,  and  Cq.  Outputs  are  printouts,  storage  on  cards 
or  files,  and  plots. 


88 


NRL-1 


THE  NRL  PARABOLIC  EQUATION  PROGRAM  PACKAGE  FOR  ACOUSTIC  PROPAGATION 


John  S.  Perkins  and  Ralph  N.  Baer 
Naval  Research  Laboratory 
Washington,  DC  20375 


We  solve  the  parabolic  equation  using  the  split-step  technique.  The 
program  makes  full  use  of  the  vectorizing  and  pipelining  capabilities  of 
the  Texas  Instruments  Advanced  Scientific  Computer  (ASC)  at  NRL.  It  is 
available  to  the  Navy  scientific  community  through  the  Navy  Laboratory 
Computer  Network  (NALCON).  In  our  implementation  we  have  included  the  op¬ 
tion  to  partially  correct  for  the  parabolic  approximations  [J.  A.  DeSanto, 
J.  S.  Perkins,  and  R.  N.  Baer,  J.  Acoust.  Soc.  Am.  64,  1664-1666  (1978)]. 

In  order  to  begin  the  range-stepping  algorithm,  two  preliminary  pro¬ 
grams  are  used  to  specify  the  environment  and  generate  an  initial  pressure 
field.  The  first  program  reads  environmental  data,  interpolates  sound- 
speed  profiles,  writes  the  bathymetry  and  a  list  of  profiles  on  a  file, 
and  draws  a  prof ile/bathyraetry  plot  and  a  sound-speed  contour.  The  second 
program  writes  a  file  containing  the  initial  pressure  field.  The  user  can 
choose  from  three  alternatives:  (1)  a  normal  mode  calculation,  (2)  a  func¬ 
tional  form  which  is  Gaussian  in  depth,  or  (3)  supply  a  complex  FORTRAN 
function  which  will  override  the  Gaussian.  A  plot  of  the  initial  field  is 
also  generated. 

In  the  main  range-stepping  program,  the  sound  speed  in  the  bottom  is 
constant  in  depth,  but  can  vary  with  range.  An  option  allows  the  bottom 
sound  speed  to  vary  so  that  a  constant  critical  angle  is  maintained.  Below 
the  water/bottom  interface,  the  square  of  the  index  of  refraction  is  given 
an  exponentially  increasing  Imaginary  part  which  effectively  produces  an 
absorbing  boundary  at  the  bottom  of  the  transform  region. 

The  direct  outputs  of  the  .sain  program  are  a  listing  and  plot  of 
transmission  loss  versus  range  for  selected  depths  (the  intensities  can  be 
smoothed  with  a  Gaussian  average  prior  to  converting  transmission 

loss),  and  a  solution  file.  Otner  plots  are  available  through  auxiliary 
plotting  programs  which  read  the  solution  file:  (1)  additional  transmis¬ 
sion  loss  versus  range  plots,  (2)  transmission  loss  versus  depth  plots, 
(3)  histograms  of  transmission  loss  values  in  a  specified  range-depth 
region,  end  (4)  gray-scale  contour  plots  of  intensity  as  a  function  of 
range  and  depth. 

It  is  also  possible  to  simulate  the  performance  of  a  tilted  array  -'f 
hydrophones  in  the  calculated  pressure  field.  This  is  done  in  an  auxtllsiy 
program  which  provides  either  gray-scale  or  Isometric  plots  of  intensity 
as  a  function  of  arrival  angle  and  range  from  the  source,  and  estimates  of 
array  signal  gain  and  3  dB  width. 
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NRL-2  AND  NRL-3 


VECTORIZED  PE  PROGRAMS 
H.  Brock 

Naval  Research  Laboratory 
Washington,  DC  20375 


The  current  PE  program  implements  the  solution  of  two  different 
parabolic  approximations.  The  "standard"  (Leontovich-Fock)  parabolic 
approximation  (NRL-2)  is  solved  using  the  Tappert-Hardin  split-step 
Fourier  algorithm.  The  second  approximation,  Tappert's  Co  independent 
equation  (NRL-3)  is  solved  with  a  finite  difference  algorithm  that 
attempts  to  match  interface  conditions  at  boundaries.  The  finite  differ¬ 
ence  algorithm  is  currently  under  development,  and  a  umber  of  problems 
have  not  been  solved.  The  code  makes  use  of  the  TI  ASC  vector  hardware 
with  explicit  in-line  generation  of  vector  instructions  and,  hence,  is  not 
readily  transportable  to  other  machines. 

The  split-step  code  is  basically  a  vectorized  version  of  the  program 
described  in  NORDA  Technical  Note  12,  The  AESD  Parabolic  Equation  Model, 
with  the  following  modifications  and  additions: 

1.  Data  structures  have  been  modified  for  more  efficient 
vector l zat ion. 

2.  The  sine  transforms  of  the  real  and  imaginary  parts  of  the  field 
are  performed  in  parallel  to  reduce  overhead  and  generate  longer  vectors* 
The  fast  Fourier  transform  algorithm  is  a  vector  version  that  eliminates 
the  bit  reversal  permutation  sequence  that  is  necessary  In  most  scalar 
implementations . 

3.  Two  source  functions  are  available — a  Gaussian  or  the  filter 
function  developed  by  SAI  (Caron,  Hanna  and  Rost). 

4.  The  pseudo  radiation  condition  boundarv  1$  Imposed  '•  y  using  a  low 
pass  digital  filter  to  attenuate  the  last  quarter  of  the  transtorm  grid. 

5.  Triangular  sector  sound  speed  interpolation  is  used  - 

b.  Two  range  step  sire  estimates  are  jutde  each  stop.  !n  physical 
space  the  range  step  size  estimate  Is  based  on  the  absolute  value  of  the 
leading  error  term  (third  order  In  range  step)  resulting  from  the  assump¬ 
tion  that  the  exponential  operators  commute.  In  Fourier  space  a  range  step 
estimate  is  derived  from  the  range  interference  length  corresponding  to 
the  Lig:.est  vertical  ^nvenunber  component  50  dii  below  the  spectral  peak. 
(The  Fourier  space  estimate  was  derived  by  H.  Garun  of  SAl.)  The  code 
selects  the  minimum  of  the  two  step  size  estimates. 
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NUSC-1 


NUSC  IFD*  PE  MODEL 

♦Implicit  Finite-Difference 

D.  Lee  and  G.  Botseas 
Naval  Underwater  Systems  Center 
New  London,  CT  06320 


A  PE  propagation  loss  model 

parabolic  wave  equation  is  solved  by  an  IFD  scheme 

The  c^de  treats  the  propagation  problem  as  an  IBVP  (Initial  Boundary 
Value  Problem) 

Sound  speed  profile  can  be  supplied  in  user's  format 

Some  bottom  boundary  conditions  are  built-in  (sloping  bottom  case) 

User's  bottom  boundary  conditions  are  acceptable 

Automatic  handling  of  horizontal  interfaces 

FORTRAN  (NUSC  VAXU  computer) 

The  workshop  test  problems  were  solved  by  this  model  using  a  Gaussian 
starter  and  the  artificial  bottom  technique. 
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NUSC-2 


THE  SNAPE  NORMAL  MODE  MODEL  FOR  THE  PARABOLIC  EQUATION 

David  H.  Wood 

Naval  Underwater  Systems  Center 
New  London,  CT  06320 


A  normal  mode  computer  model  for  predicting  how  sound  travels  in  the 
ocean  generates  solutions  of  the  Helmholtz  equation.  Such  models  can  be 
easily  modified  to  instead  solve  the  Parabolic  equation  by  changing  a  few 
lines  of  instructions.  Such  a  modified  normal  mode  model  is  useful  for 
validation  of  other  Parabolic  equation  models  based  on  other,  more  gener¬ 
ally  applicable,  numerical  techniques.  This  eliminates  two  difficulties 
usually  encountered  in  such  validation:  the  Helmholtz  equation  has  a  dif¬ 
ferent  solution  than  the  Parabolic  equation,  and  the  latter  is  sensitive 
to  a  parameter  that  does  not  even  enter  into  normal  mode  computations.  A 
normal  mode  model  can  be  further  modi  fled  tc  use  any  desired  starting 
values,  which  eliminates  another  source  of  variation  in  the  outputs  of 
Parabolic  equation  models.  Roughly  speaking,  one  merely  replaces  the 
functions  K00 )  (k  r  cos  0n)  that  express  the  range  dependence  of 
each  term  of  the  normal  mode  sum  in  terms  of  the  horizontal  wave  number 
k  cos  0n.  These  functions  are  replaced  by  the  functions  exp  -1(1/2) 
k^rCl-sin^O^ .  As  an  example  the  SNAP  model  (SACLANTCEN  Normal- 
Mode  Propagation  model)  of  Finn  B.  Jensen  and  M.  C.  Ferla  based  on  the 
work  of  A.  V.  Newman  and  F.  Ingenito  was  modified  by  Susan  M.  Bates  of 
Science  Consultants  Incorporated  of  Newport,  R.I.,  and  renamed  SNAPE. 
Numerical  results  for  some  PE  workshop  problems  have  been  generated  for 
both  Gaussian  initial  data  and  mode-sum  initial  data.  These  examples  and 
further  details  are  given  in  NUSC  TM  No.  811080. 


SACLANT 


PAREQ 

F.  B.  Jensen  and  H.  R.  Krol 
SACLANT  ASW  Research  Centre 
LaSpezia,  Italy 


PAREQ(1»2)  Is  a  parabolic  equation  model;  the  computer  program  used 
in  this  study  is  a  modified  version  of  the  one  developed  by  AESD(3), 
This  model  not  only  handles  a  variable  profile  in  depth  and  range  but  also 
allows  the  bottom  depth  and  bottom  structure  to  vary  in  range .  The  bottom 
is  characterized  by  a  compressional  speed  profile,  density  (the  density 
discontinuity  is  smoothed  using  a  hyperbolic  tangent  function  as  suggested 
in  Ref .  1)  and  attenuation  which  is  Included  by  using  a  complex  sound 
speed .  The  second  layer  of  the  bottom  has  constant  acoustic  properties* 
There  are  two  options  for  the  initial  field:  Gaussian  source  or  normal 
modes.  This  model  is  not  only  resident  on  a  UNIVAC  1106  but  also  ruA3  on 
an  HI'  21MX  computer. 

1.  Tappert,  F.  D.  The  parabolic  approximation  method.  In:  Keller,  J.  B. 
and  Papadakis,  J.  S.  eds  Wave  Propagation  and  Undorwater~Acoustics.  Lec¬ 
ture  Noces  in  Physics  70.  New  York,  Springer-Yerlag,  1977:  224-287. 

2.  Jeneen,  F.  B.  and  Krol,  H.  R.  The  use  of  the  parabolic  equation  method 
in  sound  propagation  modeling,  SACLANTCEN  SM-72,  La  Spezla,  Italy,  SAC¬ 
LANT  ASW  Research  Centre,  1975* 

3.  Brock,  H.  K.  The  AESD  parabolic  equation  nodel,  NORDA  TN-12.  NSTL 
Station,  Hiss.,  Naval  Oceen  Raserch  and  Development  Activity,  1978. 
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SAI-1 


SAI  PE 

C.  Spofford,  L.  Dozier,  H.  Garon,  and  R.  Stieglltz 
Science  Applications,  Inc* 

1710  Goodridge  Dr. 

P.  0.  Box  1303 
McLean,  VA  22102 


SAI  PE  is  based  on  the  original  AESD  PE  version  written  by  Harvey 
Brock.  Changes  in  the  original  version  implemented  in  SAI’s  PE: 

1)  A  new  source  function.  The  original  Gaussian  beam  source  function 
left  too  much  energy  traveling  at  steep  angles  and  not  enough  energy  trav¬ 
el  _  lg  at  lower  angles  at  the  source.  The  new  PE  source  function  is  an  im¬ 
pulse  response  in  physical  space  (approximately  a  sinx/x)  which  transforms 
to  a  filter  function  in  vertical  wave  number  space.  This  filter  function 
is  very  close  to  the  wave  number  spectrum  of  a  limited  aperture  point 
source  with  most  of  the  energy  evenly  distributed  over  the  low  wave 
numbers  and  a  roll-off  in  the  higher  wave  numbers  to  prevent  aliasing. 

2)  A  new  range  step  calculation.  This  feature  of  SAI's  PE  was  instal¬ 
led  to  save  computer  time.  Basically  PE  now  calculates  the  maximum  angle 
of  energy  propagation  and  then  determines  the  interference  length  (hori¬ 
zontal)  between  energy  traveling  at  this  angle  and  energy  traveling 
horizontally.  This  Interference  length  is  the  new  range  step.  This  algo¬ 
rithm  avoids  the  necessity  of  calculating  the  third  order  truncation  error 
term  (which  is  too  time  consuming  on  our  non-vec tori zing  machine)  and  pro¬ 
vides  the  same  results. 

3)  CMOD.  CMOD  is  a  transform  of  the  environment,  both  the  sound  speed 
profiles  and  the  inp"t  and  output  depths,  in  order  to  reduce  the  parabolic 
phase  velocity  error. 

4)  Bottom  options.  SAI  PE  regular  allows  the  user  to  specify  either  a 
sound  speed  profile  and  attenuation  function  in  a  sediment  layer  or  a 
transmission  loss  versus  grazing  angle  curve. 
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SAI-2 


HIGH  ANGLE  PE 

Robert  R.  Greene 
Science  Applications,  Inc* 
1710  Goodridge  Dr* 

P.  0*  Box  1303 
McLean,  VA  22102 


High  Angle  PE  is  a  full-wave  range  dependent  propagation  model 
similar  in  concept  to  PE.  It  is  designed  to  handle  the  dual  problem  of 
bottom-interaction  and  high  angles  of  propagation  up  to  about  45  degree  . 
It  is  based  upon  an  implicit  tridiagonal  range  step  and  a  cubic  spline 
discretization  in  the  vertical  coordinate. 

The  most  Important  observation  about  the  method  is  that  since  it  uses 
a  finite  difference  approach  to  the  problem  rather  than  the  split-step  al¬ 
gorithm,  a  far  better  approximation  to  high  angle  propagation  is  used  in 
the  form  of  a  "rational  parabolic"  as  opposed  to  "parabolic"  approxlma- 
atlon. 


The  method  also  allows  an  accurate  treatment  of  strong  gradients  In 
8r*md  speed  and  discontinuities  in  sound  speed  and  density.  These  are  the 
niv-in  difficulties  for  PE  in  modeling  bottom  interaction.  Strong  gradients 
are  not  a  limitation  for  a  finite  difference  technique  because  the  solu¬ 
tion  does  not  depend  on  splitting  an  operator  into  nearly  commuting  parts. 
Discontinuities  in  sound  speed  and  density  are  treated  accurately  by 
matching  cubic  splines  at  the  interfaces.  The  boundary  conditions  of  the 
wave  equation  are  the  natural  matching  conditions  for  the  splines. 

The  major  application  of  the  method  will  be  over  regions  where  bottom 
interaction  and  range  dependent  bathymetry  are  driving  the  problem.  Howev¬ 
er,  it  can  be  used  to  model  long-range  propagation.  It  should  be  particu¬ 
larly  appropriate  for  continental  slope  and  shallow  water  applications. 

The  work  la  currently  being  funded  by  ONR  Code  486. 
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IMPEDANCE  +  IPD  COMPUTER  MODEL 
John  S.  Papadakis 

Dept. .  Math.  University  of  Rhode  Island 
Kingston,  RI  02881 


Test  problems  1-4  were  run  using  an  experimental  computer  model, 
consisting  of  the  implicit  finite-difference  model  presently  under 
development  by  D.  Lee  and  G.  Botseas  NUSC,  New  London  Lab.  into  which  an 
"impedance"  subroutine  has  been  incorporated.  1:  should  be  pointed  out 
that  the  above  IMPEDANCE+IFD  computer  model  is  in  the  experimental  stage 
and  under  continuous  evolution. 

The  impedance  subroutine  calculates  the  field  at  the  advanced  range 
along  the  bottom  interface  via  a  weighted  sum  of  its  normal  derivative 
along  the  interface  from  range  zero  to  the  present  range.  This  eliminates 
the  need  of  introducing  an  artificial  bottom  in  the  IFD  model,  and  the 
field  is  advanced  only  in  the  water  column.  The  impedance  subroutine  can 
trace  the  field  along  horizontal  bottom  boundaries  as  well  as  sloping 
bottoms. 

ADVANTAGES  of  the  impedance  approach  versus  introduction  of  an 
absorbing  layer  and  a  false  bottom:  (1)  Problems  because  of  jump  discon¬ 
tinuities  along  the  bottom  interface  of  density  and  sound  speed  are  elim¬ 
inated  since  the  domain  of  interest  terminates  at  the  real  bottom.  (2)  The 
question  of  how  deep  the  false  bottom  should  be  placed  is  also  eliminated. 
(3)  The  number  of  mesh  points  in  the  z-direction  is  much  smaller  (in  many 
cases  half,  since  often,  in  the  false  bottom  approach  the  water  depth  is 
doubled).  (4)  The  computer  time  is  substantially  Improved  slp.ee,  by  (a) 
above  a  larger  range  step-size  can  be  used,  and  by  (b)  the  number  of 
z-mesh  points  Is  smaller. 

DISADVANTAGES:  (1)  If  Che  whole  weighted  sura  is  used  then  a  larger 
storage  space  is  needed.  However,  in  most  applications  because  of  the  size 
of  the  weights,  the  sum  is  truncated  and  a  constant  number  of  terms  is 
used  throughout  the  computation.  (2)  Because  of  the  use  of  a  sun,  when  the 
intensity  is  very  small  along  Che  bottom  boundary,  numerical  noise  is 
introduced.  This  can  be  seen  In  Test-Problem  1.  Actually,  the  same  disad¬ 
vantage  exists  in  the  split-step,  namely;  since  every  point  in  the  advan¬ 
ced  range  is  calculated  via  a  weighted  sum  oi  the  field  at  the  present 
range,  if  the  intensity  is  down  to  a  few  significant  figures  numerical 
noise  will  be  introduced.  This  problem  does  not  appear  in  finite  differ¬ 
ence  approach  since  a  point  at  the  advanced  range  is  given  in  terms  of  not 
all  Che  points  in  the  present  range  as  in  the  split-step  but  in  terms  of 
only  three  points.  In  other  words  the  impedance  as  well  as  the  split-step 
being  approximations  of  Integral  operators  are  global  operators  with  the 
above  disadvantage  where  the  finite  difference  being  an  approximation  of  a 
differential  operator  is  a  local  operator. 
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SACLANTCEN  NORMAL-MODE  ACOUSTIC  PROPAGATION  MODEL 


F.  B.  Jensen  and  M.  C.  Ferla 
SACLANT  ASW  Research  Centre 
La  Spezia,  Italy 


SNAP  (1)  is  a  normal-mode  model  based  on  a  program  originally  de¬ 
veloped  at  the  U.S.  Naval  Research  Laboratory  (2,  3);  this  program  solves 
the  eigenvalue  problem  by  direct  numerical  integration  of  the  depth  depen¬ 
dent  equation.  Computation  time  for  some  of  the  key  subroutines  has  been 
reduced  and  the  program  has  been  restructured  to  run  interactively  on  a 
UNIVAC  EXEC  8  system.  The  model  allows  for  slight  range  dependence  by 
employing  the  adiabatic  approximation.  SNAP  was  originally  designed  for  a 
shallow  water  environment  but  a  version  now  exists  which  handles  high- 
frequency  deep-water  situations.  Environmental  inputs  are:  arbitrary  ve¬ 
locity  profile  as  function  of  depth  (multiple  profiles  for  range-dependent 
adiabatic  computations)  in  the  water  column,  density,  attenuation  and 
compresslonal  velocity  profile  of  sediment  layer  and  density,  shear  and 
compresslonal  velocity  and  shear  and  compresslonal  attenuation  of  the 
basement.  At  present,  the  existing  output  options  are:  loss  vs.  range, 
loss  vs.  depth,  depth-averaged  loss  vs.  range,  depth-averaged  loss  vs. 
frequency,  contoured  loss  vs.  depth  and  range,  contoured  loss  vs.  fre¬ 
quency  and  range,  contoured  depth-averaged  loss  vs.  frequency  and  range, 
modal  group  velocity  vs.  frequency,  modal  phase  velocity  vs.  frequency, 
mode  function  vs.  depth,  phase  of  field  vs.  depth.  Intensity  of  field  vs. 
arrival  angle,  and  sound  speed  vs.  depth.  The  model  has  been  compared  with 
other  models  (4)  and  has  been  tested  against  various  data  sets  (5). 

1.  Jeasen,  F.  B.  and  Ferla,  M.  C.  SNAP:  the  SACLANTCEN  normal-mode 
acoustic  propagation  model,  SACLANTCEN  SM-121.  La  Spezia,  Italy,  SACLANT 
ASW  Research  Centre,  1979  (AD  A  067256). 

2.  Miller,  J.  and  Ingenlto,  F.  Normal  mode  FORTRAN  programs  for  cal¬ 
culating  sound  propagation  in  the  ocean,  NRL  Memorandum  Rpt  3071,  Washing¬ 
ton,  D.C.,  U.S.  Naval  Research  Laboratory,  1975. 

3.  Ingenlto,  F.,  Ferris,  R. ,  Kuperaan,  W.  A.  and  Wolf,  S.N.  Shallow  water 
acoustics,  summary  report  (first  phase),  NRL  Rpt  8179,  Washington  D.C., 
U.S.  Naval  Research  Laboratory,  1978. 

4.  Jensen,  P»  B.  and  Kuperaan,  W.  A.  E  * ' '  fe$*  n;:al  acoustical  modeling  at 
SACLANTCEN,  SACLANTCEN  SR-34,  U  Spezia,  Italy,  SACLANT  ASW  Research  Cen¬ 
tre,  1979  (AD  A  081  853). 

5.  Ferla,  M.  C.f  Dreinl,  G.,  Jersen,  F.  B.  and  Kuperaan,  W.  A.  Broadband 
model/data  comparisons  for  acoustic  propagation  in  coastal  waters.  Ijk  Ku- 
perman,  W.  A.  aRd  Jensen,  F.  B.  eds.  Bottom-Interacting  Ocean  Acoustics. 
New  York,  Plenum  Press,  I960. 
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PAST  FIELD  PROGRAM 


H.  W.  Kutschale 
Columbia  University 
Palisades,  NY 


FFP  (1,  2)  is  a  complete  numerical  solution  of  the  range  independent 
wave  equation  and  hence  includes  the  continuous  part  of  the  spectrum  and 
the  near  field  (distances  greater  than  a  wavelength).  The  model  used  in 
this  study  was  developed  at  Columbia  University  (2).  The  model  Includes 
propagation  of  both  compressional  and  shear  waves  and  is  therefore  also 
suitable  for  seismic  studies.  In  the  original  version,  the  depth  dependent 
part  of  the  wave  equation  was  solved  using  the  Thompson-Haskell  matrix 
method  where  the  water  column  is  divided  into  Isovelocity  layers.  An 
updated  version  (3)  allows  constant  gradient  in  (equivalent  to  l/c^ 
constant  gradient)  and  therefore  employs  Airy  functions  in  the  Thompson- 
Haskell  matrices  Instead  of  trigonometric  functions.  The  Inputs  are  the 
same  as  that  of  the  normal-mode  model  (SNAP).  There  are  two  possible  out¬ 
puts:  there  is  the  "integrand"  plot,  which  is  essentially  a  plot  of  energy 
vs.  wavenumber  and  it  therefore  has  maxima  corresponding  to  the  normal- 
mode  eigenvalues  in  the  discrete  part  of  the  spectrum.  The  maxima  in  the 
continuous  portion  of  the  spectrum  correspond  to  the  so-called  "virtual 
modes."  The  second  output  is  loss  vs.  range;  it  should  be  noted  that  the 
present  program  requires  a  complete  new  run  for  a  change  of  source  or  re¬ 
ceiver  depth.  Hence,  loss  contours  over  range  and  depth  are  not  practical 
as  compared  to  normal  mode  or  PE  methods. 

1.  DiNapoll,  F.  R.  Fast  field  program  for  multilayered  media,  NUSC  Rept. 
4103.  New  London,  Conn.,  Naval  Underwater  Systems  Center,  1971. 

2.  Kutschale,  H.  W.  Rapid  computation  by  wave  theory  of  propagation  loss 
in  the  Arctic  Ocean,  Rept.  CU-8-73.  Palisades,  N.Y.,  Columbia  University, 
1973. 

3.  Kutschale,  H.  V.  Unpublished  communication. 


98 


UNCLAS 


SECURITY  CLASSIFICATION  OF  THIS  PAGE  (TWiwi  Polo  gn(«fdj 

I  REPORT  DOCUMENTATION  PAGE 


t.  REPORT  NUMBER 


NORDA  Technical  Note  143 _ ^ 

4.  TITLE  (ond  Subtitle) 

NORDA  Parabolic  Equation  Workshop 


2.  30 VT  ACCESSION  NO. 


READ  INSTRUCTIONS 
BEFORE  COMPLETING  FORM 

3.  RECIPIENT'S  CATALOG  NUMBER 


5.  TYPE  OF  REPORT  »  PERIOD  COVERED 


Final 

6.  PERFORMING  ORG.  REPORT  NUMBER 


17.  AUTHOR/#; 

J.  A.  Davis 
D.  W.  White 
R.  Cavanagh 


»•  PERFORMING  ORGANIZATION  NAME  AND  ADDRESS 

Naval  Ocean  Research  k  Development  Activity 
Numerical  Modeling  Division,  OSTL 
NSTL  Station,  Mississippi  39529 


II.  CONTROLLING  OFFICE  NAME  AND  ADDRESS 

Same 


.  MONITORING  AGENCY  name  a  AOORESS (tldllloront  from  Controlling  O  fit  co) 


(.  CONTRACT  OR  GRANT  NUMBER/#; 


10.  PROGRAM  ELEMENT.  PROJECT.  TASK 
AREA  «  WORK  UNIT  NUMBERS 

63785N 


12.  REPORT  DATE 

September  1982 


13.  NUMBER  OF  PAGES 

102 


18.  SECURITY  CLASS,  (ol  thlo  roport) 

UNCLASSIFIED 


15a.  DECLASSIFICATION/ DOWNGRADING 
SCHEDULE 


1«.  DISTRIBUTION  STATEMENT  (ol  thlo  Roport) 

Unlimited 


17.  DISTRIBUTION  STATEMENT  (ol  tho  obotroot  on  to  tod  In  Block  30,  II  dlltoronl  tram  Roport) 


It.  KEY  WORDS  fConUnu*  on  r«vwi«  ifdt  i/ nicoiMr  end  Identity  by  block  m mber) 

Parabolic  approximation 

Workshop 

Ocean  acoustics 


0.  ABSTRACT  (Continue  on  side  If  nee eetery  md  Identify  by  biook  number) 

A  Parabolic  Equation  (PE)  workshop,  sponsored  by  the  Surveillance  Environmen¬ 
tal  Acoustic  Support  (SEAS)  Project  and  hosted  by  NORDA's  Numerical  Modeling 
Division,  was  held  at  NSTL  from  31  March  to  3  April  1981.  The  purpose  of  the 
workshop  was  to  provide  a  forum  for  those  active  in  theoretical  and  applied  PE 
development  and  to  compare  computer  results  for  a  set  of  ocean  acoustic  prob¬ 
lems.  Fifteen  formal  presentations  were  given,  and  thirteen  different  PE  models 
were  exercised  against  four  test  cases.  The  test  cases  ranged  from  one  that  any 


FORM 
I  JAN  73 


EDITION  OF  I  NOV  68  13  OBSOLETE 
S/N  01 OMF-01 4-4601 


UNCLASSIFIED 

SECURITY  CLASSIFICATION  OF  THIS  PAGE  /WA#fl  Doto  Bntorod) 


tKCUMTY  CLASSIFICATION  OF  THIS  FAOE  (Wlim  Of  Kntf«0 


PE  should  be  able  to  handle  to  one  that  is  a  challenge  for  every  PE .  This  report 
describes  the  results  of  that  workshop  and  includes  the  theoretical  bases  for 
the  models,  the  model  descriptions,  results  of  the  test  cases,  and  abstracts  of 
the  formal  presentations. 


UNCLASSIFIED 


ucurnrv  classification  of  yhu  i 


DEPARTMENT  OF  THE  NAVY 

OFFICE  OF  NAVAL  RESEARCH 
875  NORTH  RANDOLPH  STREET 
SUITE  1425 

ARLINGTON  VA  22203-1995 


IN  REPLY  REFER  TO: 


5510/1 

Ser  321OA/01 1/06 
31  Jan  06 


MEMORANDUM  FOR  DISTRIBUTION  LIST 

Subj :  DECLASSIFICATION  OF  LONG  RANGE  ACOUSTIC  PROPAGATION  PROJECT 
(LRAPP)  DOCUMENTS 

Ref:  (a)  SECNAVINST  55 10.36 

Enel:  (1 )  List  of  DECLASSIFIED  LRAPP  Documents 

1 .  In  accordance  with  reference  (a),  a  declassification  review  has  been  conducted  on  a 
number  of  classified  LRAPP  documents. 

2.  The  LRAPP  documents  listed  in  enclosure  (1)  have  been  downgraded  to 
UNCLASSIFIED  and  have  been  approved  for  public  release.  These  documents  should 
be  remarked  as  follows: 

Classification  changed  to  UNCLASSIFIED  by  authority  of  the  Chief  of  Naval 
Operations  (N772)  letter  N772A/6U875630,  20  January  2006. 

DISTRIBUTION  STATEMENT  A:  Approved  for  Public  Release;  Distribution  is 
unlimited. 

3.  Questions  may  be  directed  to  the  undersigned  on  (703)  696-4619,  DSN  426-4619. 


f 

BRIAN  LINK 
By  direction 


Subj :  DECLASSIFICATION  OF  LONG  RANGE  ACOUSTIC  PROPAGATION  PROJECT 
(LRAPP)  DOCUMENTS 

DISTRIBUTION  LIST: 

NAVOCEANO  (Code  N121LC  -  Jaime  Ratliff) 

NRL  Washington  (Code  5596.3  -  Mary  Templeman) 

PEO  LMW  Det  San  Diego  (PMS  181) 

DTIC-OCQ  (Larry  Downing) 

ARL,  U  of  Texas 

Blue  Sea  Corporation  (Dr.Roy  Gaul) 

ONR  32B  (CAPT  Paul  Stewart) 

ONR  3210A  (Dr.  Ellen  Livingston) 

APL,  U  of  Washington 

APL,  Johns  Hopkins  University 

ARL,  Penn  State  University 

MPL  of  Scripps  Institution  of  Oceanography 

WHOI 

NAVSEA 

NAVAIR 

NUWC 

SAIC 


Declassified  LRAPP  Documents 


